Geostatistically estimation and mapping of forest stock in a natural unmanaged forest in the Caspian region of Iran

Document Type: Research Paper

Authors

1 Research Institute of Forests and Rangelands Tehran

2 Science and Research Branch Tehran

3 University of Tehran

Abstract

Estimation and mapping of forest resources are preconditions for management, planning and research. In this study, we applied kriging interpolation of geostatistics for estimation and mapping of forest stock at-tributes in a natural, uneven-aged, unmanaged forest in the Caspian region of northern Iran. The site of the study has an area of 516 ha and an elevation that ranges from 1100 to 1450 m a.s.l. Field sampling was per-formed on a 75m × 200m systematic grid using 309 geo-referenced circular sample plots of 1000 m2 area. Experimental variograms were calculated and plotted for basal area (BA), volume (V) and stem density (N). Whereas the calculated variograms of BA and V exhibited spatial auto-correlation only after data stratification based on diameter size classes and tree species, the variogram of stem density displayed a moderate spatial structure that was fitted by a spherical model. Stem density was estimated by ordinary block kriging and the accuracy of estimation was validated by cross-validation result. We conclude that geostatistical approaches have the potential to more accurately capture and describe the spatial variability of forest stock, and thus reduce the uncertainty in estimates of stem density as well as produce more accurate stem density maps of forests in comparison with the spatially uninformed classic method. Geostatistical methods provide a very suitable tool to derive more accurate estimates of growing stock, particularly in structurally complex, unmanaged, uneven-aged forest such as this one from the Caspian region of northern Iran.

Keywords


Geostatistically estimation and mapping of forest stock in a natural unmanaged forest in the Caspian region of Iran

 R. Akhavan1*, H. Kia-Daliri2, V. Etemad3

 1- Forest Research Division, Research Institute of Forests and Rangelands, Tehran, Iran

2- Faculty of natural Resources, Islamic Azad University, Science and Research Branch, Tehran, Iran

3- Dept. of Forestry, University of Tehran, Karaj, Iran

 * Corresponding author’s E-mail: akhavan@rifr-ac.ir

(Received: April. 07.2014   Accepted: July. 30.2014)

ABSTRACT

Estimation and mapping of forest resources are preconditions for management, planning and research. In this study, we applied kriging interpolation of geostatistics for estimation and mapping of forest stock at-tributes in a natural, uneven-aged, unmanaged forest in the Caspian region of northern Iran. The site of the study has an area of 516 ha and an elevation that ranges from 1100 to 1450 m a.s.l. Field sampling was per-formed on a 75m × 200m systematic grid using 309 geo-referenced circular sample plots of 1000 m2 area. Experimental variograms were calculated and plotted for basal area (BA), volume (V) and stem density (N). Whereas the calculated variograms of BA and V exhibited spatial auto-correlation only after data stratification based on diameter size classes and tree species, the variogram of stem density displayed a moderate spatial structure that was fitted by a spherical model. Stem density was estimated by ordinary block kriging and the accuracy of estimation was validated by cross-validation result. We conclude that geostatistical approaches have the potential to more accurately capture and describe the spatial variability of forest stock, and thus reduce the uncertainty in estimates of stem density as well as produce more accurate stem density maps of forests in comparison with the spatially uninformed classic method. Geostatistical methods provide a very suitable tool to derive more accurate estimates of growing stock, particularly in structurally complex, unmanaged, uneven-aged forest such as this one from the Caspian region of northern Iran.

 

Keywords: Geostatistics, mapping, forest stock, unmanaged forest, Caspian region, Iran.  

 


INTRODUCTION

The implementation of close-to-nature silvi-culture that aimed at conserving and restor-ing unique, important and ecologically and economically endangered (due to human dis-turbances and illegal exploitation)  Caspian forests of Iran (Abdollahpour & Assadi Atui, 2005) that resemble natural forests, requires a thorough understanding of ecological pro-cesses that have created the structures and composition of these natural archetypes. For-tunately, natural mixed uneven-aged hard-wood deciduous forests that encompass all three major development stages (i.e., the de-cay, initial, and optimal stages; Leibundgut, 1993, Korpel, 1995) still exist along the southern part of the Caspian Sea and the north-facing aspects of Elborz Mountain in northern Iran (Sagheb-Talebi et al., 2004). These Caspian forests are characterized by a high plant species diversity of more than 80 woody species (Sagheb-Talebi et al., 2004) and afford a unique opportunity to study the spatial structures found in intact natural old-growth forests. Although proper maps of spatial patterns of trees are a prerequisite for forest planning and management, planners rarely have wall-to-wall coverage of spatially explicit tree stem maps that depict the spatial distribution of forest stock attributes (e.g., stem density, basal area, and volume) over large areas. Instead, forest stock attributes are typically inventoried using a sampling process that lacks spatial resolution of stock at-tributes between sample points. Further, non-spatial independence suggests that many statistical tools and inferences may not be appropriate and a spatial extrapolation of results to any given point in the whole sample area may not be warranted. Thus, an un-biased interpolation method has to be employed for mapping forest stock attributes after sampling that permits forest managers to more closely delineate the major development stages based on natural spatial structures of tree species and enable researchers to more realistically model stand structures and dynamics. Geostatistics, which is concerned with the detection, modeling and estimation of spatial dependence of continuously distributed or regionalized variables (Isaak & Srivastava, 1989, Goovaerts, 1997), is a useful tool to describe the spatial structure and pro-vide geo-referenced estimates of forest attributes. Although the spatial distribution of trees in a particular stand represents a point pat-tern of discrete objects (Dale, 2000), stand at-tributes such as basal area, volume stock and density can be thought to be directly influenced by different spatially continuous varia-bles such as solar radiation, soil characteris-tics and water nutrient availability, thus al-lowing to be considered as spatially continuous (Kint et al., 2003). Geostatistics, thus, pro-vides a natural framework for estimation techniques in forest inventory sampling (Mandallaz, 1991) and has frequently been used to estimate and map forest resources based on forest-scale surveys (Samra et al., 1989, Biondi et al., 1994, Gunnarsson et al., 1998, Tuominen et al., 2003, Montes et al., 2005, Freeman & Moisen, 2007, Tröltzsch et al., 2009, Brus et al., 2012), enabling the production of geo-referenced data maps for basal area, density, standing volume and tree species groups at scales where these variables often reveal spatial auto-correlation. Following two studies in which geostatistics was applied to estimate and map forest stock at-tributes in a natural managed forest(Akhavan et al., 2010) and a plantation forest (Akhavan & Kia-Daliri, 2010) in the Caspian region of Iran, this study aimed at extending this approach to a natural, unmanaged forest in this region of Iran. This study specifically aimed at (1) investigating the spatial variability and spatial structure of forest stock attributes (stem density, basal area, and volume) using variogram analysis in geostatistics, (2) estimating and mapping forest stock attributes using the geostatistical kriging interpolation method, and (3) examining whether kriging improves the estimation accuracy of these forest attributes in comparison with a spatially uninformed classic approach.

 

MATERIALS AND METHODS

Study Area

Caspian forests along the southern part of the Caspian Sea, span the elevational gradient of north-facing aspects of Elborz Mountain in northern Iran, occur in five main vegetation types (i.e., Querco-Buxetum, Querco-Carpinetum, Parrotio-Carpinetum, Fagetum hyrcanum and Carpinetum orientale) and are remnant of the Tertiary era that extend for 800 km in an east to west direction (Sagheb-Talebi et al., 2004). The study area (516 ha) is located in the northern part of the fourth district (Chelir) of the educational and research forest station of Tehran University (Kheiroud) at 51˚ 40' E lon-gitude and 36˚ 32' N latitude (Fig. 1). Eleva-tion varies from 1100 m to 1450 m above sea level and slopes range from 5% to 85%. The mean annual temperature, precipitation and relative humidity are 15.3 ˚C, 1458 mm, and 83%, respectively. The climate is cold and wet in the winter and temperate in the summer without any dry season. The growing season is 270 days per year. The forest under investi-gation is a typical naturally seeded, uneven-aged, mixed hardwood old-growth stand in the Fagetum hyrcanum association that experi-enced no silvicultural activity and has no known management history. The inventoried forest area is a mixture of broad-leaf decidu-ous tree species. The dominant tree species are beech (Fagus orientalis Lipsky) and horn-beam (Carpinus betulus L.), alongside with maple (Acer velutinum Boiss. and A. cappadoci-cum Gled.), alder (Alnus subcordata C.A.M.)and oak (Quercus castaneifolia C.A.M.). Large topographic variation contributes to the heterogeneous nature of the stand (Fig. 1).

 

 

 

Fig. 1. Study area and sampling grid (75 m × 200 m).

 


Data collection

A network with 75 m (N-S) × 200 m (W-E) systematic rectangular grid was used for sapling. At each grid point, we established a circular sample plot of 1000 m2 surface area where we recorded the UTM coordinates of the plot center (Fig.1) and the diameter at breast height (d.b.h. in 1.3 m above the ground) of each tree with a d.b.h. that exceeded 7.5 cm. Further, the height of the closest tree to its respective plot center point and the largest tree in each plot were measured for volume estimation. The inventory was accomplished in the summer of 2011. Stand-level attributes of interest, i.e., basal area (BA), volume stock (V) and stem density (N) were computed for all trees as well as for trees in different size classes (i.e., small size [d.b.h. ≤ 32.5 cm, S], medium size [32.5< d.b.h. ≤ 52.5, M], large size [52.5 < d.b.h. ≤72.5 , L], and extra large size [d.b.h. > 72.5 cm, EL]; after Sagheb-Talebi et al., 2005), and separately for the two most dominant tree species (i.e., beech and hornbeam) that accounted for about 75% of the total species frequency in the study area.

 

Geostatistical approach for spatial auto-correlation analyses

To identify and describe the spatial depend-ency (i.e., spatial auto-correlation) of BA, V and N of all trees and the diameter size clas-ses of the two dominant tree species, we used variogram (semi-variance) analysis as our geostatistical approach. Three parameters are commonly used to describe and model the behavior of variogram: range, sill and nugget effect. The range is the distance where the spatial correlation disappears and the variogram levels off, the sill corresponds to the height of the variogram after leveling off, and the nugget effect is represented by the intercept of variogram on the ordinate axis. The ratio of the nugget effect to the sill is known as the relative nugget effect. This is a measure of the percentage of the variability in the data from sources other than spatial auto-correlation. A low relative nugget effect (≤ 25%) is a sign of strong spatial auto-correlation where the application of geostatistical techniques is particularly beneficial. Nugget effects between 25% and 75% indicate moderate spatial auto-correlation and high relative nugget effect (≥75%) can indicate either weak spatial autocorrelation in the population or spatial patterns at scales smaller than the sampling distance (Cambardella et al., 1994, Ganawa et al., 2003, Freeman & Moisen, 2007). In the current study, all the experimental variograms obtained were modeled using either a pure nugget effect or an additional spherical model. The selection criterion used was the minimal residual sum of squares. The spherical model is given by (Webster & Oliver, 2000): 

2000):   for 

 for                               (1)

Where h, c0, c, and a represent a particular lag vector, nugget effect, structural variance, and range, respectively. After normalization of the data, only omnidirectional (isotropic) variograms were calculated, because experi-mental variogram surfaces showed no vario-gram anisotropy. By common convention, the analysis was restricted to distances of half of the study area dimension (i.e., 2200 m). We further used the Kriging interpolation method as our geostatistical approach for possible interpolating values among sample plots and mapping the distribution of BA, V and N. Kriging computes surfaces of the best linear unbiased estimation of regionalized variables at un-sampled points based on the spatial structure defined by the experimental semi-variogram. Ordinary kriging (the most common type of kriging in practice, particularly in environmental sciences) of the regionalized variable   at point     is given by (Webster & Oliver, 2000): 

                                         (2)

Where,  is the weight associated with the value of at the sampled point  with the nonbiased condition:

                                                      (3) Kriging may be used for estimation at a single point (point kriging) or over an area (block kriging). In this study, since the mean is as-sumed to be stationary and unknown and because no large-scale trend was observed, ordinary block kriging without trend was used. A 32 m × 32 m mesh (approximately the same area as a sample plot to emphasize the local variation around the sampling plots) was used to discretize the study area for block kriging interpolation. The estimates were ob-tained using the nearest 16 data plots within the maximum effective range of variograms that corresponded to the scale of auto-correlation. A cross-validation (leave-one-out) was performed for selected datasets to evaluate the kriging results. Cross-validation was evaluated by calculating the Mean Bias Error (MBE), which should ideally be equal zero, because kriging is unbiased (Webster & Oli-ver, 2000):

                             

 (4)

                                   (5)

The accuracy of kriging was measured using the Root Mean Square Error (RMSE):

              (6)

                              (7)

 Where  is the estimated value of re-gionalized variable at the location of ; N is

The number of sample plots and  is the mean value of measured samples of interest-ed attribute. The software package used for geostatistical analysis was GS+ version 9 (Gamma Design Software, LLC, Plain Well, MI).

 

Point pattern analyses

To further understand the spatial dependen-cies of BA, V and N, we computed several indices for quadrat counts (e.g., variance-to-mean ratio, Morisita’s index of dispersion, and Morisita’s standardized index of disper-sion) to analyze separately the spatial point patterns for all trees, trees in different diame-ter size classes, and for beech and hornbeam. The variance-to-mean ratio, attributed to Fisher et al. (1922) is one of the oldest and simplest measures of dispersion. The ratio        ( ) usually called the index of disper-sion (I) (Bailey & Gatrell, 1995) and is based on the observation in a random pattern, de-scribed by the Poisson distribution, the vari-ance equals the mean, so I = 1 for a random pattern. Ratios larger than 1 indicate clump-ing, while smaller ratios indicate regular or uniform pattern. The other frequently used quadrat-based dispersion index is the Stand-ardized Morisita index of dispersion. Smith-

 

Gill (1975) set out to improve Morisita's index (Morisita, 1962) just described by putting it on an absolute scale from -1 to +1. To calculate this index, Morisita's index of dispersion (Id) was first calculated, along with two critical values of the uniform index (Mu) and the clumped index (Mc):

                                  (8)

                                (9)

                               (10)

 

Where n is the sample size, x is the number of individuals,  and  are the values

 

of chi-squared with (n-1) degrees of freedom that have 2.5% or 97.5% of the area to the

right. Morisita's Id is 1 for a random distribution, >1 for a clumped distribution, and

a regular/uniform distribution. The Standardized Morisita index (Ip) is then calculated by one of the four following formulae:

 

When

      (11)

When         (12)

When        (13)

When 

  (14)

The Ip ranges from -1.0 to +1.0, with 95% con-fidence limits at +0.5 and -0.5. Random pat-terns give an Ip of zero, clumped patterns give an index value of above zero and uniform patterns below zero.

 

RESULTS

Sampling

A total of 309 plots were sampled in the study area. Normalization tests showed that the data were not normally distributed, conse-quently, the data were transformed using the square root and log transforms to more closely approximate normal distributions. Table 1 shows the summary statistics of the sample plots. The distribution of tree sizes in the study area reveals an uneven-aged stand diameter distribution with the majority of trees in the small and medium tree size classes (Table 2) and a few very large trees up to 300 cm d.b.h. (Fig. 2).

 

Geostatistical analyses

Because we did not find any variogram aniso-tropies, we fitted only omni-directional vario-grams using a spherical model to which a nugget effects was added (Fig. 3, Table 3). No spatial auto-correlation was revealed in the experimental variograms for BA and V. Con-sequently, these attributes represented a pure nugget effect. However, N showed a moderate level of spatial auto-correlation (50%).

 

 

 

 

 

Table 1. Summary statistics of 309 sample plots

Attribute

Mean

Min

Max

SD

CV

%

Skewness

Kurtosis

BA (m2 ha-1)

35.8

8.6

99.6

11.5

32.2

1.34

5.13

V (m3 ha-1)

459.7

95.1

1452.4

191.7

41.7

1.37

4.15

N (n ha-1)

279

20

1140

166.2

59.6

1.61

3.86

SD: Standard Deviation; CV: Coefficient of Variation

 

 

Table 2. Frequency and proportion of trees in the four diameter classes in the study area

Relative frequency

Absolute frequency

Diameter class

64.6%

5561

S

17.4%

1493

M

9.8%

850

L

8.2%

708

EL

100%

8612

Total

S, Small size (d.b.h. ≤ 32.5 cm); M, Medium size (32.5

L, Large size (52.5 72.5 cm)

 

 

Table 3. Characteristics of fitted variograms for the three attributes of basal area (BA), volume (V), and stem density (N)

Attribute

Lag

(m)

Fitted

model

Nugget

effect

Sill

 

Range

(m)

R2

%

SpD

 %

BA (m2 ha-1)

130

Pure nugget effect

0.86

0.89

-

7

96.6

V (m3 ha-1)

130

Pure nugget effect

0.11

0.12

-

37

91.6

N (n ha-1)

130

Spherical

0.11

0.22

533

93

50

SpD (Spatial Dependence) = (SpD ≥75%, Weak; 25%< SpD <75%, Moderate; SpD ≤ 25%, Strong)

 

 

Fig. 2.Diameter distribution of trees in the study area.

 

 

 

Because BA and V did not show any spatial auto-correlation, kriging interpolation was used only for N (Table 4). A comparison of tables 1 and 4 reveals that the estimated mean stem densities are not very different. However,

 

 

compared to the spatially uninformed classic method, a variance reduction in the estimate of N of approximately 72% was achieved using kriging interpolation. Figure 4 shows a kriging map and an error map for stem density over the study area.

 

 

Fig. 3. Graphs of isotropic variograms and fitted models for the three attributes of basal area (BA), volume (V), and stem density (N) in the study area. Solid lines represent the models, filled circles represent point-estimates. The minimum number of pairs for each lag. distance is 274

 

 

 

 

Table 4. Summary statistics of kriging results for forest stem density (n ha-1)

Mean

Min

Max

SD

CV%

260

67

597.32

82.5

31.7%

.

 

Evaluating the kriging results of N with the Mean Bias Error and Root Mean Square Error and their respective relative values (Table 5)

revealed a low relative Mean Bias Error of less than 10% (MBEr = 6.9%). A cross-validation graph of stem density confirms the accuracy of our estimation (Fig. 5).

 

Variography after data classification based on tree size classes and species

Although BA and V did not show any spatial variability when considering all trees, the re-sults were quite different when we computed the experimental variograms by tree size clas-

 

ses and species, individually. When consider-ing trees of different diameter size classes, the spatial variability of BA (Fig. 6) and V (Fig. 7) clearly weakens from the small diameter size class towards the extra large class, with an approximately moderate spatial dependence observed in the small and medium size classes and only a pure nugget effect in the large and extra large classes. Similarly, when examining the spatial variability of BA by tree species, both beech and horn-beam individually exhibited moderate spatial auto-correlation (Fig. 8) that was completely masked when examining the spatial variability of all species (Fig. 6)

 

 

Table 5. Validation results of kriging interpolation for stem density (n ha-1)

MBE

RMSE

MBEr

RMSEr

19.4

147.7

6.9%

52.9%

MBE: Mean Bias Error; RMSE: Root Mean Square Error; MBEr: relative MBE; RMSEr: relative RMSE

 

 

 

 

Fig. 4. Kriging map (a) and error (standard deviation) map (b) of stem density in the study area

 

 

 

 

Fig. 5. Cross-validation graph for stem density (N). The solid line represents the regression and the dotted line is the 45- degree line; the best unbiased estimation occurs when the solid line coincides with the dotted line


 

 

Point pattern analyses

Point pattern analyses revealed that the spatial distributions of all trees, the two dominant tree

 

 

species and all diameter size classes except the extra large class were clumped in the study area (Table 6).

 

 

 

Fig. 6. Graphs of isotropic variograms and fitted models for basal area (BA) based on the four diameter size classes in the study area. Solid lines represent the models, filled circles represent point-estimates. Minimum number of pairs for each lag distance is 258. Note: S= Small size, M= Medium size, L= Large size, and EL= Extra Large size.

 

Fig. 7. Graphs of isotropic variograms and fitted models for volume (V) based on the four diameter size classes in the study area. The abbreviations and explanations are the same as for Fig. 6.

 

Table 6. Spatial patterns for all trees, size classes and the two dominant tree species in the study area

 

Index

 

 

 

Id

Ip

Spatial patterns

All trees

9.91

1.32

0.500

Clumped*

S

15.91

1.83

0.501

Clumped*

M

2.6

1.33

0.500

Clumped*

L

1.58

1.21

0.500

Clumped*

EL

0.958

0.98

-0.139

Uniform*

Beech

6.47

2.01

0.501

Clumped*

Hornbeam

14.86

1.87

0.501

Clumped*

variance-to-mean ratio; Id, the Morisita’s index; Ip, the Standardized Morisita’s index; *, significant at p

 

 

Fig. 8. Graphs of isotropic variograms and fitted models for basal area (BA) for the two most abundant tree species in the study area. Solid lines represent the models, filled circles represent point-estimates. Minimum number of pairs for each lag distance is 274.

,


DISCUSSION

The spatial structure of the three forest attrib-utes (i.e., basal area, volume and stem density), expressed by their spatial autocorrelation, differed in this unmanaged natural uneven-aged deciduous forest. Whereas stem density behaved as a regionalized variable and exhibited a moderate spatial structure (Table 3), overall stand basal area and volume did not show any spatial auto-correlation and exhibited a pure nugget effect (Fig. 3). The reason for these apparent differences can be found in the spatial distribution of trees of different diameter size classes. While smaller trees exhibited a clumped spatial distribution, larger trees increasingly tended toward a random and regular spatial distribution. This development toward spatial randomness or regularity with increasing tree size has been demonstrated in many forest types (e.g., Szwagrzyk & Czerwczak, 1993, Zenner & Peck, 2009) and is the main reason why the kriging interpolation method embedded in a geostatistical approach was able to estimate stem density more accurately to reduce the variance for stem density estimation by approximately 70%, and to produce a smaller coefficient of variation with acceptable estimation accuracy (MBEr ≈ 7%; Table 5) in comparison with the spatially uninformed classic approach (Table 1). Because smaller trees were more numerous than larger trees in this old-growth forest (Fig. 2) and each individual tree contributes equally to stem density, the moderate spatial structure of the stem density of the stand is largely driven by the

 

clumped spatial pattern of the smaller trees (Table 6). This is in stark contrast to the observed lack of spatial structure of the stand basal area and volume (Fig. 3) that are principally driven by the spatial randomness of the larger trees that, despite their low density, account for a large proportion of the growing stock. Thus, the spatial distribution of larger trees has more influence on the spatial structure of the growing stock (basal area and volume) in forests than that of small and medium-sized trees that typically exhibit a moderate spatial structure for basal area, and volume (Figs. 6 and 7). This lack of spatial structure of basal area and volume appears to be independent of changes induced by management, as was shown for a stand in the nearby Namkhane district (Fig. 1) that was subject to a management plan that prescribed a close-to-nature silviculture (Akhavan et al., 2010) and resulted in a diameter distribution that was similar to that of the unmanaged Chelir forest investigated in this study. Because larger trees were retained in the Namkhane forest following two harvest entries (Fig. 9), no spatial structures were detected for stand basal area and volume in Namkhane as well. The spatial structure of stem density differed between the unmanaged and man-aged stands, however, and was most likely influenced by management interventions such as road construction, plantation, and harvesting that may have masked a clear signal of the spatial structure of stem density in the Namkhane district (Akhavan et al., 2010).

 

 

 

Fig. 9. Diameter distribution of trees in the Namkhane district; after Akhavan et al. (2010).

 

 

The spatial structures of stem density and basal area in the current study are quite different from those reported by Akhavan and Kia-Daliri (2010), who applied the same geo-statistical approach in an eighteen year-old maple (Acer velutinum Boiss.) plantation in the Caspian region of Iran. In that study, stem density did not exhibit any spatial structure whereas basal area did. This was likely due to the regular spacing of the trees and the homogeneity of tree sizes in the plantation that was still in the early stage of stand development. Neither windstorms and droughts nor intra-specific competition had been strong enough to remove stems and break up the initial regular spatial distributions (initial planting space of 3 m × 3 m) after 18 years. Consequently, the regularity of the spatial distributions of the plants remained intact and no spatial auto-correlation was found for stem density. In contrast, the homogeneity of tree sizes in the maple plantation had a low coefficient of variation for d.b.h (35%), which was less than half of the coefficient of variation for d.b.h. observed in this study (79%), was sufficient to induce spatial auto-correlation for stand basal area. However, it is expected that as stem numbers decline until harvesting time (at approximately 80 years), stem density will exhibit spatial auto-correlation and the spatial structure of basal area will disappear as the homogeneity of tree sizes decreases. The low, medium and high density areas that are typically associated with mature, middle-aged, and young stands, respectively, have become clearly visible in the kriged stem density maps (Fig. 4a). These kriged stem density maps, thus, permit an indirect estimation of forest stock for any point in the area. For example, the point I on the digital map of figure 4a has a stem densi-ty of 170 n ha-1. Based on the diameter size classification in table 2, there are 64.6% (110 trees), 17.4% (30 trees), 9.8% (16 trees), and 8.2% (14 trees) of trees in the small, medium, large and extra large size classes at the point I, respectively. Therefore, if we use the mid-point value of each size class (namely, 20, 42.5 and 62.5 cm for the first three size classes, respectively; for the largest size class this depends on the maximum diameter size in the studied area), we can indirectly estimate the growing stock at point I. To obtain a more precise estimate, more and narrower size classes in the upper range of the size distribution would of course be beneficial due to their relatively larger impact on forest volume stock. Nonetheless, the advantage of a geostatistical approach that enables the precise quantification of the statistical error in a kriged error map (Fig. 4b) is that it identifies areas with higher errors that can then be covered with extra sample plots to reduce this estimation error. Each tree species has a specific physiological age and longevity that is reflected in a specific auto-correlation diagram (i.e., variogram) that may be masked when investigating the overall stand structure. While we failed to detect any spatial au-to-correlation for basal area when analyzing all tree species that were present at the same time (Fig. 3), we did detect a moderate spatial structure for basal area for the two dominant species (i.e., beech and hornbeam with clumped spatial distribution; Table 6), individually (Fig. 8). Hence, it appears that the spatial structures of growing stock attributes in mixed stands, which are often uneven-aged in this region, may be weaker than those observed in pure and even-aged stands, particularly if not analyzed separately by species. For this reason, Akhavan and Kia- Daliri (2010) may have been able to detect a spatial structure for basal area in their plantation forest. In general, we hypothesize that if the diameter distribution observed in a natural forest has a narrower range (e.g., no larger trees present) and has no missing size classes (i.e., is without any interruption), the growing stock appears to have a stronger spatial structure with a higher auto-correlation, but this tentative hypothesis needs to be put to more rigorous testing. Geostatistical approaches were instrumental for the identification of the spatial structures of stem density in this natural unmanaged forest. In close-to-nature forestry, which is increasingly being applied throughout the world, kriged stem density maps can provide greater assistance for forest management planning and for the identification of forest development stages (Akhavan et al., 2012) than traditional tree size class maps. Whereas the identification of develop-ment stages is straightforward and fairly easy in monospecific and even-aged forests, the detection of these stages is often difficult and imprecise in mixed and uneven-aged forests due to different habitat requirements of tree species. Here, kriged stem density maps could make an effective contribution toward a more accurate identification and delineation of these stages. Spatial interpolation of forest tree attributes using geostatistical approaches (kriging) can provide useful outputs at the regional scale as well, facilitating rapid spatial analyses that enable the user to test the influence of different factors on forest stock attributes (Merganič et al., 2004). These maps can be very useful for forest owners and man-agers as a suitable tool to gain more insight into the spatial distribution of the growing stock. This information is also beneficial as a guide map for writing forest management plans, i.e., to identify harvesting areas and to plan road networks, inter-planting, and silvi-cultural interventions (i.e., thinning and light thinning) that are based on forest density and stock distribution, to ensure that forestry is done more sustainably.

 

ACKNOWLEDGMENTS

This research was financially supported by Iran National Science Foundation (INSF) with grant No. 88001186 and partially supported by the Research Institute of Forests and Rangelands (RIFR), Iran. The authors thank M. Hasani for assistance with fieldwork and data collection and Kh. Mirakhorlou for assistance with the preparation of maps.

Abdollahpour, M. & Assadi Atui, A.R. (2005). Entwicklung der geregelten Forstwirtschaft in den Kaspischen Wäldern. In: Nosrati K., Marvie Mohadjer, R., Bode, W. and Knapp, H.D. (eds). Schutz der Biologischen Vielfalt und integriertes Management der Kaspischen Wälder (Nordiran), pp. 99-121.

Akhavan, R. & Kia-Daliri, H. (2010). Spatial variability and estimation of tree attributes using geostatistics in a plantation forest in the Caspian region of Iran. Caspian Journal of Environmental Sciences, 8: 163-172.

Akhavan, R., Zahedi Amiri, Gh. & Zobeiri, M. (2010). Spatial variability of forest growing stock using geostatistics in the Caspian region of Iran. Caspian Journal of Environmental Sciences, 8: 43-53.

Akhavan, R., Sagheb Talebi, Kh., Zenner, E.K. & Safavimanesh, F. (2012). Spatial patterns in different forest development stages of an intact old-growth Oriental beech forest in the Caspian region of Iran. European Journal of Forest Research, 131: 1355-1366.

Bailey, T.C. & Gatrell, A.C. (1995). Interactive spatial data analysis. Prentice Hall, New York, 413 pp.

Biondi, F., Myers, D.E. & Avery, C.C. (1994). Geostatistically modeling stem size and increment in an old-growth forest. Canadian Journal of Forest Research, 24: 1354-1368.

 Brus, D.J., Hengeveld, G.M., Walvoort, D.J.J., Goedhart, P.W., Heidema, A.H., Nabuurs, G.J. & Gunia, K. (2012). Statistical mapping of tree species over Europe. European Journal of Forest Research, 131: 145-157.

Cambardella, C.A., Moorman, T.B., Parkin, T.B., Karlen, D.L., Turco, R.F. & Konopka, A.E. (1994). Field scale variability of soil properties in central Iowa soil. Soil Science Society of America Journal, 58, 1501-1511.

Cressie, N.A.C. (1993). Statistics for spatial data. John Willy and Sons, Inc., New York, 900 p.

Dale, M.R.T. (2000). Spatial pattern analysis in plant ecology. Cambridge University press, Cambridge, United Kingdom. 326 pp.

Fisher, R.A., Thornton, H.G. & Mackenzie, W.A., 1922. The accuracy of the planting method of estimating the density of bacterial population. Annals of Applied Biology, 9: 325-359.

Freeman, E.A. & Moisen, G.G. (2007). Evaluating kriging as a tool to improve moderate resolution maps of forest biomass. Environmental Monitoring and Assessment, 128: 395-410.

Ganawa, E.S.M. Soom ,M.A.M. Musa, M.H., Mohammad Sharif, A.R. & Wayayok, A. (2003). Spatial variability of total nitrogen and available phosphorus of large rice field in Sawah Sepadan Malaysia. Science Asia, 29: 7-12.

Goovaerts, P. (1997). Geostatistics for natural resources evaluation. Oxford University Press, New York 483 p.

Gunnarsson F., Holm S.,Holmgren P. and Thuresson T. (1998). On the potential of kriging for forest management planning. Scandinavian Journal of Forest Research, 13: 237- 245.

Isaaks, E.H. & Srivastava, R.M. (1989). An introduction to applied geostatistics. Oxford University Press, New York, 561 pp.

Kint, V., Meirvenne, M.V., Nachtergale, L., Geudens, G. & Lust, N. (2003). Spatial methods for quantifying forest stand structure development: a comparison between nearest neighbour indices and variogram analysis. Forest Science, 49: 36-49.

Korpel, S. (1995). Die Urwälder der Westkarpaten. Gustav Fischer, Berlin, Germany, 310 p.

Leibundgut, H. (1993). Europäische Urwälder. Verlag Paul Haupt Bern and Stuttgart, Germany, 260 p.

Mandallaz, D. (1991). A unified approach to sampling theory for forest inventory based on infinite population model. Ph.D. thesis, academic press, ETH Zürich, Switzerland, chair of forest inventory and planning. 257 p.

Mandallaz, D. (1993). Geostatistical methods for double sampling schemes: application to combined forest inventory. Technical report, ETH Zürich, chair of forest inventory and planning, 133 p.

Merganič, J., Quednau, H.D. & Šmelko, Š. (2004). Relations between selected geomorphology features and tree species diversity of forest ecosystems and interpolation on a regional level. European Journal of Forest Research, 123: 75- 85.

Morisita, M. (1962). "Iδ-Index, A Measure of Dispersion of Individuals". Researches on Population Ecology, 4: 1–7.

Montes, F., Hernandez, M.J. & Canellas, I. (2005). A geostatistical approach to cork production sampling in Quercus suber forests. Canadian Journal of Forest Research, 35: 2787-2796.

Sagheb-Talebi, Kh., Sajedi, T. & Yazdian, F. (2004). Forests of Iran. Research Institute of Forests and Rangelands press, Technical publication, No. 339, 28 p.

Samra, J.S., Gill, H.S. & Bhatia, V.K. (1989). Spatial stochastic modeling of growth and forest resource evaluation. Forest Science, 35: 663-676.

Smith-Gill, S.J. (1975). Cytophysiological basis of disruptive pigmentary patterns in the leopard frog, Rana pipiens. II. Wild type and mutant cell specific patterns. Journal of Morphology, 146: 35-54.

Szwagrzyk, J. & Czerwczak, M. (1993). Spatial pattern of trees in natural forests of East-Central Europe. Journal of Vegetation Science, 4: 469-476.

Tröltzsch, K., van Brusselen, J. & Schuch, A. (2009). Spatial occurrence of major tree species groups in Europe derived from multiple data sources. Forest Ecology and Management, 257: 294-302.

Tuominen, S., Fish, S. & Poso, S. (2003). Combining remote sensing, data from earlier inventories, and geostatistical interpolation in multi-source forest inventory. Canadian Journal of Forest Research, 33: 624- 634.

Webster, R. and Oliver, M.A. (2000). Geostatistics for environmental scientists. Wiley press. 271 p.

Zenner, E.K. & Peck, J.E. (2009). Characterizing structural conditions in mature managed red pine: spatial dependency of metrics and adequacy of plot size. Forest Ecology and Management, 257: 311–320.