Spatial analysis of carbon storage density of mid-subtropical forests using geostatistics:a case study in Jiangle County,southeast China

Zhuo Lin?Lin Chao?Chengzhen Wu?Wei Hong?Tao Hong?Xisheng Hu

1 Introduction

The global carbon cycle is a natural process that occurs in Earth’satmosphere and biosphere and involvesthe exchange of matter and energy(Falkowski et al.2000).Forests are considered to be the main carbon sinks in the terrestrial ecosystem and play an important role in the carbon cycle(Fang et al.2001).Forest carbon storage is the amount of carbon that is retained by vegetation after assimilation through photosynthesis and autotrophic respiration.In recent years,forest carbon storage has received considerable interest by researchers,mainly due to its capacity of offsetting CO2emissions to the atmosphere.Many activities affect forest carbon storage,such as forest growth,expansion,harvesting,and burning(Granier et al.2000).However,little is known about the precise impact of these factors or how to quantitatively evaluate the distribution of forest carbon storage(Fang et al.2001;Pan et al.2011;Saatchi et al.2011).Previous studies have focused on calculating forest carbon storage by traditional methods,which has not included an effective examination of spatial characteristics.Moreover,the relationship between the determinants and carbon storage density(CSD)has yet to be fully described.The distribution of forest carbon storage is not uniform across geographical regions due to signi ficant differences in species composition and richness,topographic features,and vegetative cover.Forest carbon storage may have high spatial heterogeneity and spatial autocorrelation,with continuous variation on the spatiotemporal scale.

The mid-subtropical forest is one of the largest components of the subtropical forest in China.Because subtropicalforests,whose complexity rankssecond to tropical forests,have abundant natural resources and biodiversity,much mid-subtropicalforesthas been exploited and destroyed.Most of the existing forests are secondary or plantation forests.Although a number of studies have reported CSD in tropical and temperate forests,and the carbon storage of mid-subtropical forest has received attention in recent years in China,CSD in mid-subtropical remains unclear(Zhang et al.2007).Geostatistics has been considered an effective approach to facilitate the quantification of spatial variation(Matheron 1963)and is widely used to analyze spatial heterogeneity and to assess spatial autocorrelation distribution(Montes et al.2005).Spatial analysis has been used in forestry,where some variables exhibit the characteristics of spatial arrangements(Hirata et al.2009;Herna′ndez-Stefanoni et al.2011;Lamsal et al.2012).A variety of geostatistical techniques,such as semi-variogram,Kriging,and autocorrelation analysis,have been applied to explore spatial characteristics in forestry in recent years.Keane et al.(2012)quanti fi ed spatial variability of fuel-component biomass using spatial variograms and found that each component has its own inherent scale in low-elevation northern Rocky Mountain forest.Mora and Beer(2013)used geostatistics to uncover spatial relationships between root length density of Coffea arabica and soil nutrientrelated factors at the plot scale.However,few studies have considered the spatial distribution of forest carbon storage(Liu et al.2014).This is especially true for midsubtropical forest,despite its important role in in fl uencing carbon emissions and sinks.

Thus,the main object of our study was to quantitatively evaluate the geo-spatial distribution of mid-subtropical forest carbon storage.We used the typical mid-subtropical forests of Jiangle County as a study case.Geostatistics was used to depict geographic variations of CSD,a key indicator of forest carbon sequestration.We first used inventory data from sample plotsto modelthe spatial heterogeneity of CSD via a semi-variogram analysis at the village level.Next,we used the parameters of the optimal theoretical semi-variogram model in a Kriging interpolation to map the spatial distribution of CSD for the entire county.Subsequently,global and local spatial autocorrelation analyses determined the characteristics of the CSD spatial pattern.Then,we employed regression analysis methods to explore the factors affecting CSD in Jiangle County.Finally,based on the results,forest management measures were proposed to increase carbon storage in the mid-subtropical forest of Jiangle County.

2 Methods

2.1 Study area

This study was conducted in Jiangle County(117°05′E to 117°40′E and 26°26′N to 27°04′N)of Sanming City in Fujian Province,southeast China(Fig.1).Fujian Province has the highest forest coverage in China(approximately 65%)and the most complete mid-subtropical forest ecosystem.Jiangle County is one of China’s most important forestry areas and has a total forest area of over 189,000 hectares(approximately 98.4% ofthe mountain lands).The forest coverage rate is 85.2%,ranking the county first in Fujian Province.The forests of this county have typical characteristics of mid-subtropical forest,including a variety of important coniferous timber;mixed forests;and complete,natural,subtropical,evergreen broad-leaved forests in the Longqi Mountain Nature Reserve,whose forest coverage reaches 98%.The terrain of the county is rugged;geographical types are primarily hills and low mountains,with an average slope of 27°.The elevation is mainly below 800 m.The study area is under a subtropical monsoon climate with characteristics of an oceanic climate and continental climate.The annual average temperature is 18.7°C,the annual average precipitation is 1698.2 mm,the annual average sunshine time is 1730 h,and the average frost-free period is 295 days.

Fig.1 The study area.a The location of Fujian Province in China;b the location of Jiangle County in Fujian Province;c map of Jiangle County showing its 13 towns

The forest resources of Jiangle County can be divided into four types based on elevation:I.Evergreen broad-leaf forest zones.These regions are typical vegetation zones that range from about 200–1000 m elevation.Major arbor tree species include:Castanopsis canesii(Castanopisis carlesii(Hemsl.)Hayata),Chinkapin(Castanopsis sclerophylla(Lindl.)Schott),Lithocarpus(Lithocarpus skanianus(Dunn)Rehd.),Quercus glauca(Cyclobalanopsis glauca(Thunb.)Oerst.),and Red nandu(Machilus thunbergii Sieb.et Zucc.);II.Mixed broadleaf-conifer forest zones.These regions are relatively small transition zones between evergreen broad-leaf and coniferous forest zones.Major arbor tree species include:Masson pine(Pinus massoniana Lamb),Schima superba(Schima superba Gardner&Champ.),and Mocketprivetlike oak(Quercus phillyreoides A.Gray);III.Coniferous forest zones.These regions host major commercial species,including Chinese fi r(Cunninghamia lanceolata(Lamb.)Hook.),Masson pine(Pinus massoniana Lamb),and Cryptomeria fortunei(Cryptomeria japonica var.sinensis Miq.);.Moso bamboo forest zones.These regions have Moso bamboo (Phyllostachys pubescens)along with broad-lea for conifer species and cover a total area of about 300 km2.Major arbor tree species include:Chinese fi r(Cunninghamia lanceolata(Lamb.)Hook.),Masson pine(Pinus massoniana Lamb),Paulownia (Paulownia fortunei(Seem.)Hemsl),Wild jujube(Ziziphus jujuba Mill.var.spinosa(Bunge)Hu ex H.F.Chow),Birch(Betula luminifera H.Winkl.),Sassafras(Sassafras tzumu(Hemsl.)Hemsl.),Chinese tulip (Liriodendron chinense (Hemsl.)Sarg),Camphor tree(Cinnamomum camphora(L.)J.Presl),Sweet-scented osmanthus(Osmanthus fragrans(Thunb.)Lour.),Manglietia yuyuanensis(Manglietia yuyuanensis Y.W.Law),Mieheliamaudiae(Michelia maudiae Dunn,Alniphyllum fortunei(Hemsl.)Makino),Liquidambar(Liquidambar formosana Hance),Choerospondias axillaris(Choerospondias axillaris(Roxb.)B.L.Burtt& A.W.Hill),Phoebe bournei(Phoebe bournei(Hemsl.)Yen C.Yang), Elaeocarpus sylvestris (Elaeocarpus sylvestris(Lour.)Poir.),and Ormosia henryi(Ormosia henryi Prain).

2.2 Data collection

The stand and tree data used in this study were collected from the Forest Resource Inventory of Fujian in 2012.The sample plots were distributed across all villages of Jiangle County and the total area of all survey plots reached approximately 39,000 hectares,or almost one- fi fth of thetotal forested area in Jiangle County.1The authority responsible for all plots is the Forestry Bureau of Jiangle County.The studies in these fields did not involve endangered or protected species.All activities for research were permitted by the responsible authority.Simple statistics of the research data are listed in Table 1.The locations of the 124 villages are shown in Fig.2.We considered the village as the basic analysis unit.We used the longitude and latitude of the village centroid as the basic unit coordinate.Each plot was 20 m×30 m and every tree of each plot was measured.For each plot,the following data were collected:tree species composition,diameter at breast height(DBH)(cm)of individual trees(only measured DBH≥3 cm),and individual height of all living trees(m).Other variables were computed based on the survey,including vegetation coverage(%),average age of stand(years),average DBH of stand(cm),average height of stand(m),number of trees per hectare,and volume of living trees(m3/ha).Topographic descriptors and geographic location were also recorded,and the data were calibrated using the digital elevation map(DEM)of Jiangle County provided by the Geospatial Data Cloud from the Computer Network Information Centre of the Chinese Academy of Sciences.The DEM dataset includes elevation(m),slope gradient(°),slope aspect,slope position,soil type,latitude,and longitude.

Table 1 Descriptive statistics of research data used in this study

Fig.2 Distribution of geometric centers of the 124 villages

2.3 Carbon storage density calculation

Previous studies show that the biomass expansion factor(BEF)can be used to convert measured timber volume to tonnage of carbon(Liski et al.2006;Neilson et al.2007).Therefore,we first used the volume-biomass conversion models to calculate biomass(Fang et al.2001).We then used the ratio of root and stem recommended by IPCC(Eggelston et al.2006)to get biomass of roots,since the model established by Fang et al.(2001)only includes aboveground biomass.After obtaining biomass(including aboveground and roots)for each plot,we used 0.5 as the average carbon content rate(Eggelston et al.2006).The sum of every tree species was the total carbon storage(Mg)for a given sample plot,and the carbon storage per hectare,namely,CSD(Mg/ha),was calculated as the total carbon storage of the sample plot divided by the plot area.According to survey plot data,the average CSD was 43.19 Mg/ha and standard deviation was 12.06 Mg/ha at the village level in Jiangle County.The present study only considered carbon storage of arbor species and did not include carbon storage from shrubs,herbs,soil,and litter,i.e.,the CSD in this paper only represents that of living trees,including the stems,branches,foliage and roots.Different designations(such as plot,village,town,and county)refer to areas administered by different scopes.

2.4 Geostatistics

Geostatistics is a branch of applied statistics that focuses on studying variables that are distributed continuously in space(Rossi et al.1992).Two important spatial analysis techniques are semi-variogram,which is a measure of the spatial variability of a regionalized variable,and Kriging,which provides estimates for an unrecorded location(Krige 1951;Isaaks and Srivastava 1989;Zawadzki et al.2005).

The semi-variogram is the most common measure for characterizing spatial continuity as a statistical model of the spatial dependence structure(Isaaks and Srivastava 1989).The semi-variance function used to compute a semivariogram is thus estimated at each lag distance and direction.The details of the formula and calculation processing are presented elsewhere(Fortin et al.1989;Isaaks and Srivastava 1989;Webster and Oliver 2007).In this case,the parameters of the semi-variogram were calculated using DPS 7.05(Tang and Zhang 2013)and GS+version 7.0.

The Kriging approach,originated by Krige(1951)and developed by Matheron(1963),takes into account both the distance and the degree of variation between known data points;the distribution of regional variables can be mapped through spatial interpolation based on the structure de fi ned by the fi tted theoretical semi-variogram model(Webster and Oliver 2007).In this study,we used Ordinary Kriging to produce an interpolated distribution map of CSD of Jiangle County via ArcGis 10.2.Usually,two error statistics are used to compare Kriging results with the observed values:the Kriged average error(KAE)for unbiasedness and the Kriged reduced mean square error(KRMSE)for coherence(Sarangi et al.2005).Detailed formulas and more information have been previously published(Joumel and Huijbregts 1978;Isaaks and Srivastava 1989).

Three indicators,including Global Moran’s I,Anselin Local Moran’s I,and Local Getis-Ord Gi*,were used to explore spatial autocorrelation.Global Moran’s I was introduced by Moran(1950)and examines whether spatial correlation exists over an entire region.Given a set of variables and an associated attribute,Global Moran’s I evaluates whether the pattern expressed is clustered,dispersed,or random.The value of Moran’sIranges approximately from+1 to-1,representing positive to negative spatial autocorrelation,with 0 representing the absence of spatial autocorrelation.Moran’s I statistic is a global statistic and cannot identify statistically significant spatial clusters and outliers.Ansel in Local Moran’s I can be used to explore a certain spatial regularity.The values of Anselin Local Moran’s I can indicate clusters and outliers that represent positive and negative spatial autocorrelation;the variables with no spatial autocorrelation correspond to‘‘not significant.’’Another local autocorrelation indicator is Local Getis-Ord Gi*,which considers positive spatial autocorrelation and enables differentiation between clusters of similar values that are high or low relative to the mean,thus aiding in the detection of unusual events.The detailed formulas and more information about spatial autocorrelation analysis are available for further review(Getis and Ord 1992;Anselin 1995;Ord and Getis 1995,2001;Anselin et al.2006;Nelson and Boots 2008).

2.5 Determinants regression analysis

Ordinary least squares(OLS)linear regression can be used to generate predictions or to model a dependent variable in terms of its relationships to a set of explanatory variables.Meanwhile,geographically weighted regression(GWR)can provide tests for the existence of spatial interactions between variables across locations and reveal spatial variations in empirical relationships between variables that would otherwise be ignored in OLS analysis.Therefore,these two regression models were used to evaluate the in fluence of CSD and some determinants,including local stand condition indices and topographic features,on tree and stand growth.The detailed formulas and more information about OLS and GWR can be found in the work of Brunsdon et al.(1998,2002).

In this study,a Gaussian function was determined to assign the spatial weight matrix of GWR through models constructed in ArcGIS 10.2.

3 Results

The application of geostatistics requires that the intrinsic hypothesis for a regionalized variable be met and that the raw data follow a normal distribution(Zawadzki et al.2005).CSD at the village level was not normally distributed;therefore,we processed the raw CSD data using a logarithmic transformation to explore the geographic variations.

Table 2 The parameters of spatial heterogeneity of carbon storage density for different models

3.1 Spatial heterogeneity of carbon storage density

Table 2 presents the results of four models of CSD in Jiangle County.The parameter C0/（C0+C）for each model indicates that non-structural or random factors affected CSD at the village level,but that structural factors led sections for spatial heterogeneity in the study area.Three statistical cross-validation indices were used to compare the goodness of the theoretical semi-variogram models,including the coefficient of determination(R2),the normalized mean(NM)and the normalized root mean square error(NRMSE).As a practical rule,the higher the value of R2,the better the fi tting result.In addition,results are acceptable when the absolute value of NM is close to 0 and NRMSE is close to 1.Thus,except for the linear model,the models satisfied the hypothesis of spatial heterogeneity according to the evaluation indices.Among the three satisfactory models,the Gaussian model had the highest R2(0.779)and best cross-validation results.Therefore,it was optimal to use the Gaussian model to describe the spatial distribution features of CSD in Jiangle County.The output of the Gaussian model(Table 2)showed that random factors accounted for only 25.14%of the spatial heterogeneity of CSD,whereas spatial autocorrelation accounted for approximately 74.86%.Further,the Gaussian model suggested that the observations were fully independent when the observation distance exceeded 2147 m(A0).The CSD distribution in Jiangle County had an obvious spatial autocorrelation pattern,which was affected by structural factors on the mesoscale(approximate range from 2000 m to 3000 m).Taken together,the Gaussian model accurately represented the spatial heterogeneity of CSD in Jiangle County,i.e.,the Gaussian model could be considered as the best- fitting model to obtain the CSD distribution map via Kriging.

3.2 Spatial distribution of carbon storage density

Figure 3 shows the detailed geographic distribution of CSD in Jiangle County using Ordinary Kriging based on the Gaussian model as the theoretical semi-variogram.The values of KAE and KRMSE of the cross-validation results were-0.099 and 0.886,respectively.As a practical rule,the KAE value should be close to 0 to be acceptable,whereas KRMSE should be in the range0.746 to 1.254 in this case)(Sarangi et al.2005).These guidelines indicate that the Kriging results had good predictability and consistency between the prediction errors and the standard deviation of the observed values,i.e.,the CSD distribution map of Jiangle County by Kriging was reliable and applicable.

The CSD of Jiangle County was not evenly distributed but varied greatly between regions,as shown in Fig.3a.The highest values of CSD were in the east(approximately 46 to 61 Mg/ha).In contrast,the regions of the north and northwest had relatively low CSD(approximately 34 to 39 Mg/ha).The other areas had moderate CSD(approximately 39 to 46 Mg/ha).Figure 3b shows the average Kriging interpolation values of CSD for Jiangle County’s 13 towns obtained using ArcGIS.High-CSD(＞45 Mg/ha)towns were distributed in the east and southeast of the county,including Guyong,Gaotan,Shuinan,Moyuan,and Nankou Towns.Among these,Moyuan Town had the highest CSD (55.25 Mg/ha).Trend surface analysis showed that CSD tended to increase gradually from west to east.From north to south,CSD initially increased and then declined,but the range of variation was minor in this direction.

3.3 Spatial autocorrelation of carbon storage density

In the analysis of spatial autocorrelation,the value of the neighborhood size calculated by inverse distance squared was 13651 m.The value of Global Moran’s I statistic was 0.130,and the score was 2.077(statistically significant at the 5%level),indicating that the spatial distribution pattern of CSD had a tendency towards geographical clustering of similar villages in Jiangle County.According to Ansel in Local Moran’s I(statistically significant at the 5%level),shown in Fig.4a,eight villages had positive spatial autocorrelation,including seven High–High(HH)clusters and one Low–Low(LL)cluster.Areas of negative spatial autocorrelation,including three Low–High(LH)outliersand one High-Low(HL)outlier,were also detected in the study area.Most of the HH clusters were located around Moyuan;in addition,some areas of Gaotan and Anren belonged to HH clusters.LL clusters were only detected around Dayuan.Outlier areas were mostly located in the north or northwest of Moyuan,which were indicated as HH clusters.Similarly,the results of Getis-ord G*in Fig.4b show that cold-spots were distributed mostly in the northern region of Jiangle County and were mainly located in the towns of Dayuan and Wan’an.On the basis of different significance levels,the cold-spot areas of the villages were classified into the following two groups:eight villages at the 5%significance level and nine villages at the 10%significance level.(No villages were detected at the 1%significance level.)Regarding CSD hot-spots,the results from Fig.4b were nearly the same as those detected by Ansel in Local Moran’s I in Fig.4a.The numbers of hot spot villages detected at different significance levels were five(1%),eight(5%),and five(10%),all in the eastern regions of Jiangle County.The majority of the high-CSD clusters were in Moyuan.By analyzing the results of Ansel in Local Moran’s I and Getis-ord G*,we conclude that these two statistical in dices accurately reflected the local spatial autocorrelation of CSD in Jiangle County.

Fig.3 Geographic distribution of CSD in Jiangle County by using Ordinary Kriging.a CSD distribution;b the average CSD of different towns in Jiangle County

3.4 coefficients of determinants

OLS and GWR models were established respectively to evaluate the effects of different determinants on CSD at the village level.In the models,CSD observations were the dependent variable.Because CSD(Mg/ha)was determined by many factors,a diagnostic analysis was required to determine the fitting of the independent variables that satisfied the spatial modelling assumptions.A number of stand and tree condition in dices and topographic features were tested by correlation analysis and stepwise regression.Five variables were selected as predictors based on significance testing for modelling:average DBH(cm),average tree height of stand(m),stand density(trees/ha),elevation(m),and slope(°).As OLS is a global regression approach,coefficients ofdeterminants were constant:0.065626(DBH),0.514574(height),0.001945(density),-0.002031(elevation),and 0.151959(slope).coefficients of the GWR model were considered reliable after testing multicollinearity by comparing condition numbers(smaller than 30).In modelling spatial data,analysis of the spatial autocorrelation in the regression model residual is required,and the model regression residual should be spatially random.Global Moran’s I and Z scores of the estimation residuals of OLS and GWR were used to detect these spatial characteristics.The Z score(3.31)of model residuals by OLS was larger than 2.58 at the 0.01 significance level.This indicates that coefficients of the OLS model might not be accurate because the existing significant spatial autocorrelation would cause a serious violation of the independence assumption of the OLS model and lead to inefficient estimation(Zhang et al.2005).However,the spatial autocorrelation of the GWR model residuals was not significant(Zscore of 1.55wass mall er than1.65at the0.1 significance level),indicating GWR successfully reduced the spatial autocorrelation in the model residuals and would not generate under-predictions or over-predictions(Zhen et al.2013).Thus,it was explicitly better to use GWR to explore relationships between CSD and determinants through comparative analysis of various aspects.

Fig.4 Results of local spatial autocorrelation analysis.a Cluster map of CSD detected by Anselin Local Moran’s I(statistically signi fi cant at the 5%level);b cold and hot spot distribution of CSD detected by Getis-Ord G*

Fig.5 Contour maps of coefficients of fi ve determinants estimated by GWR.a Contour map of DBH;b contour map of height;c contour map of density;d contour map of elevation;e contour map of slope.The categories of each coefficient were classified using the Quartile Method

Because GWR produced one set of model coefficients for the regression model at each location,the local spatial variation of parameters was fully incorporated in to the modelling process(Brunsdon et al.1998,2002).The spatial variation of each determinant coefficient of the GWR model is illustrated in Fig.5.This figure reveals that relationships between CSD and the five determinants varied from region to region across the study area—the larger the value of the model coefficient,the greater the influence on CSD.In general,the coefficients of DBH,height,and density were all positive over the study area,while coefficients of the elevation and slope ranged from negative to positive depending on the region(Fig.5).Although coefficients of the three stand condition indices(DBH,height,and density)were positive,the impacts of these determinants were different.The tendencies of DBH and density were similar.For example,the coefficient values of these two determinants revealed that in the southern regions of Jiangle County(covering Wanquan,Huangtan,Nankou,Moyuan,and Bailian),the influences on CSD by DBH and density were relatively higher than in other regions,as shown in Fig.5a,c.In contrast,greater influences of tree height appeared in most areas around the towns in the northeast,including Anren,Wan’an,and Gaotan(Fig.5b).The impacts on different regions were also distinct for the two topographic variables,elevation and slope.The relationship between CSD and elevation had a negative correlation in the southern part,including Wanquan,Huangtan,and Bailian(Fig.5c).Negative impacts of slope on CSD were discovered in the middle region of Jiangle County,mostly in Guyong,Shuinan,and Moyuan as well as some regions of Dayuan and Anren(Fig.5e).These results demonstrate that the magnitudes and signs of these coefficients varied by geographic region,indicating that the influences of the five determinants on CSD depended to a great extent on forest conditions and topographic features.

4 Discussion

4.1 Spatial characteristics of forest carbon storage

The distribution of forest carbon storage is heterogeneous across geographical regions(Sales et al.2007;Liu et al.2014).Thus,it is crucial and necessary to consider spatial characteristics of forest carbon storage for management planning and policy making.More advanced spatial analysis methods have been developed and applied to enhance the description of spatial distributions,such as geostatistics and GWR(Zhang et al.2009).

Geostatistics is based on the regionalized variables theory:distribution in space and spatial autocorrelation.However,some researchers recognize that geostatistics should not always be employed to analyze forest or vegetation distribution because the results might not be accurate and useful,even when accounting for spatial autocorrelation or spatial heterogeneity.Gunnarsson et al.(1998)posited that hardwood volume distribution(approximately 314 m2per plot)does not have obvious spatial autocorrelation in Sweden.Similarly,Akhavan and Kia-Daliri(2010)indicated that Kriging has no potential for simulating natural forest stock in the Caspian region of Iran.However,in the present study,these two fi ndings were contradicted:the results of the spatial heterogeneity analysis indicate that there was a relatively strong spatial autocorrelation of CSD,distributed in the mesoscale range(approximately 2000–3000 m)(Table 2),and that global and local spatial autocorrelation existed(Fig.4).The cause of this contradiction might be that spatial analysis should be based on lag distance(variation range)(Fortin et al.1989);if the spatial distance of a study area is much larger than the variation range,spatial heterogeneity and spatial autocorrelation gradually disappear with increasing spatial distance.Therefore,some studies based on the large-or global-scale could not detect spatial variations using geostatistics;determining the appropriate study scale is important for spatial analysis by geostatistics.In this case,using the village as the unit,geostatistics enabled the capture of spatial characteristics and exploration of the spatial pattern for mid-subtropical forest CSD distribution at the county level.

Determinants potentially affecting CSD included DBH,height,density,elevation,and slope.This study explored the relationships between CSD and determinants of stand and topography by using OLS and GWR.We found that some stand determinants had positive effects on CSD in Jiangle County.DBH and height had more impact on CSD than did density according to coefficients of the OLS model.This suggests that a larger volume,which usually is related to DBH and tree height,would sequester more carbon in the forests.A higher stand density could improve forest carbon storage too,but this determinant reflected tree number and had relatively little in fl uence.This finding for mid-subtropical forests is consistent with the results of other studies,such as Du et al.(2010)and Liu et al.(2014).Elevation and slope were the most important topographic factors.However,the coefficient of the OLS model for elevation was negative(-0.002031),but the absolute value of this coefficient was not very large,indicating that elevation had only a slightly negative effect on carbon storage.In contrast,the coefficient of the slope was positive(0.151959),suggesting that steeper slopes supported more carbon.

Since the OLS model violated the assumption of independently distributed errors without considering spatial autocorrelation(Legendre 1993),these constant coeffi cients did notaccurately describe the effects of determinants on forest carbon for different local areas,where forest conditions and topographic features vary.coefficients of the GWR model revealed more detailed and accurate information in the study area.The effects of these determinants on CSD varied from area to area(Fig.5).In general,DBH,height,and density had a positive correlation with CSD over the study area,but their magnitudes varied substantially by location.The values of the model coefficients of elevation and slope ranged from negative to positive,depending on the sample location.In Jiangle County,there are many planted forests in the young and mid-age stages.Some ecological processes of these forests,such as seed dispersal,seedling establishment,and survival,are in fl uenced by changes in stand conditions,including light and soil properties.Moreover,many commercial timber forests exist in some regions,particularly the north and west of the county.After harvesting,forest gap can create a variety of favorable microenvironments for regeneration.It may lead to clumped distribution of juvenile forests and negative association with mature forests(Kuuluvainen et al.1998;Are′valo and Ferna′ndez-Palacios 2003).In addition,different tree species and forest succession stages determine different spatial characteristics of CSD in mid-subtropical forests.Shade-tolerant,late-succession tree species can grow under less-clumped spatial distributions in the absence of significant exogenous disturbances(McDonald et al.2003),such as Quercus glauca forest in southwest Jiangle County(Fig.4).The presence of these different stand conditions and topographic positions tend to result in different relationships between CSD and the determinants in different regions.Thus,the localized(GWR)model coefficients could provide detailed information about the effects of determinant geographical variation on CSD.

4.2 Measures for improving forest carbon storage

Spatial analysis of CSD distribution revealed that there were geographic variations for the mid-subtropical forest of Jiangle County.Forest management measures are proposed below to improve the carbon storage of mid-subtropical forests.

The first measure is increasing the forest area.There is abundant shrub wood and open forest land in some areas of Jiangle County,leading to low CSD,such as in the northern region(Fig.3).These lands could be transformed into forestlands by scientific silviculture techniques.Moreover,to improve carbon storage,more water conservation forests should be established on steep slope land with poor site conditions.Improving the qualities of the current forests is another important measure to increase carbon storage.The CSD of mixed broad leaf-conifer forests is higher than that of pure forest(Zhang et al.2007).Therefore,afforestation measures could be implemented to optimize the structure of the tree species to increase the mixed forest areas in regions with low CSD,such as closing hillsides to facilitate afforestation and cultivating high-quality native tree species.Finally,it is necessary to make reasonable forest management decisions and develop harvesting strategies.The eastern regions of Jiangle County with high CSD(Fig.3)are currently concentrated zones for commercial forests.In these regions,a large proportion of forests are middle-aged and mature stands of timber forest,such as Chinese fir and Mass on pine.If these forests were harvested,forest carbon storage would be reduced.A rational harvesting strategy would reduce carbon storage loss.The forests of the southwest are mainly subtropical evergreen broad-leaved forests in the Longqi Mountain Nature Reserve.The broad-leaved young forest area is large,leading to a relatively lower CSD(Fig.3).Protecting natural evergreen broad-leaved forests through management would allow young stands to increase carbon storage.

5 Conclusion

This research focused on utilizing geostatistics to examine CSD spatial distribution of the mid-subtropical forest in Jiangle County.Using the semi-variogram approach to explore spatial heterogeneity,we found that the observations of CSD are dependent on the range of the lag distance of less than the mesoscale(approximately 2000 m to 3000 m).In addition,the results show that the Gaussian model had the best fitting performance.Based on this model,we mapped the overall CSD distribution of Jiangle County using Kriging analysis.The map illustrates increasing CSD from west to east of the study area.Three statistical in dices were used to analyze the spatial autocorrelation of CSD.The global indicator(Global Moran’s I)revealed that the spatial pattern of CSD was clustered(at the 5%level of significance)at the village level.Next,by calculating the local indicators(Anselin Local Moran’s I and Getis-Ord G*),we obtained more precise information on the spatial associations among the individual locations(clusters or outliers)and hot(or cold)spots for CSD in the study area.In addition,the GWR model hinted at different impacts on different areas for each of the major determinants,including DBH,height,density,elevation,and slope.Generally,DBH,height,and density had positive correlation with CSD in Jiangle County,but with a range of magnitudes.Finally,we proposed some forest management measures to improve the carbon storage of the mid-subtropical forest of Jiangle County,including increasing forest area,improving the qualities of the current forests,and implementing reasonable forest management decisions and harvesting strategies.The findings of this study provide useful information about CSD distribution for the quantification of carbon storage in the mid-subtropical forest in Jiangle County.The method of geostatistics used in this case could be used to estimate and analyze forest carbon storage for other similar regions.

AcknowledgementsThis study was jointly supported by Science and Technology Major Project of the Hall of Science and Technology of Fujian(No.2012NZ0001)and the Project of National Natural Science Fund of China(No.30671664).

Akhavan R,Kia-Daliri H(2010)Spatial variability and estimation of tree attributes in a plantation forest in the Caspian region of Iran using geostatistical analysis.Casp J Environ Sci 8:163–172

Anselin L(1995)Local indicators of spatial association—LISA.Geogr Anal 27(2):93–115. doi:10.1111/j.1538-4632.1995.tb00338.x

Anselin L,Syabri I,Kho Y(2006)GeoDa:an introduction to spatial data analysis.Geogr Anal 38(1):5–22.doi:10.1111/j.0016-7363.2005.00671.x

Are′valo JR,Ferna′ndez-Palacios JM(2003)Spatial patterns of trees and juveniles in a laurel forest of Tenerife,Canary Islands.Plant Ecol 165(1):1–10.doi:10.1023/A:1021490715660

Brunsdon C,Fotheringham AS,Charlton M(1998)Geographically weighted regression.J R Stat Soc D-Stat 47(3):431–443.doi:10.1111/1467-9884.00145

Brunsdon C,Fotheringham AS,Charlton M(2002)Geographically weighted summary statistics—a framework for localised exploratory data analysis. Comput Environ Urban 26(6):501–524.doi:10.1016/S0198-9715(01)00009-6

Du H,Zhou G,Fan Wet al(2010)Spatial heterogeneity and carbon contribution of aboveground biomass of mo so bamboo by using geostatistical theory.Plant Ecol 207(1):131–139.doi:10.1007/s11258-009-9659-3

Eggelston HS,Buendia L,MiwaK et al(2006)Guidelines for national greenhouse gas inventories.IPCC National Greenhouse Gas Inventories Programme,Hayama

Falkowski P,Scholes RJ,Boyle EEA et al(2000)The global carbon cycle:a test of our knowledge of earth as a system.Science 290(5490):291–296.doi:10.1126/science.290.5490.291

Fang J,Chen A,Peng C et al(2001)Changes in forest biomass carbon storage in China between 1949 and 1998. Science 292(5525):2320–2322.doi:10.1126/science.1058629

Fortin MJ,Drapeau P,Legendre P(1989)Spatial autocorrelation and sampling design in plant ecology.Vegetatio 83(1–2):209–222.doi:10.1007/BF00031693

Getis A,Ord JK(1992)The analysis of spatial association by use of distance statistics.Geogr Anal 24(3):189–206.doi:10.1111/j.1538-4632.1992.tb00261.x

Granier A,Ceschia E,Damesin C et al(2000)The carbon balance of a young beech forest.Funct Ecol 14(3):312–325.doi:10.1046/j.1365-2435.2000.00434.x

Gunnarsson F,Holm S,Holmgren P et al(1998)On the potential of Kriging for forest management planning.Scand J For Res 13(1–4):237–245.doi:10.1080/02827589809382981

Herna′ndez-Stefanoni JL,Dupuy JM,Tun-Dzul F et al(2011)In fl uence of landscape structure and stand age on species density and biomass of a tropical dry forest across spatial scales.Landsc Ecol 26(3):355–370.doi:10.1007/s10980-010-9561-3

Hirata A,Kamijo T,Saito S(2009)Host trait preferences and distribution of vascular epiphytes in a warm-temperate forest.Plant Ecol 201(1):247–254.doi:10.1007/s11258-008-9519-6

Isaaks EH,Srivastava RM(1989)Applied geostatistics,vol 2.Oxford University Press,New York

Joumel AG,Huijbregts CJ(1978)Mining geostatistics.Academic Press,New York

Keane RE,Gray K,Bacciu V et al(2012)Spatial scaling of wildland fuels for six forest and rangeland ecosystems of the northern Rocky Mountains,USA.Landsc Ecol 27(8):1213–1234.doi:10.1007/s10980-012-9773-9

Krige DG(1951)A statistical approach to some basic mine valuation problems on the Witwatersrand.J Chem Metall Min Soc S Afr 52:119–139

Kuuluvainen T,Ja¨rvinen E,Hokkanen TJ et al(1998)Structural heterogeneity and spatial autocorrelation in a natural mature Pinus sylvestris dominated forest.Ecography 21(2):159–174.doi:10.1111/j.1600-0587.1998.tb00670.x

Lamsal S,Rizzo DM,Meentemeyer RK(2012)Spatial variation and prediction of forest biomass in a heterogeneous landscape.J For Res 23(1):13–22.doi:10.1007/s11676-012-0228-6

Legendre P(1993)Spatial autocorrelation:trouble or new paradigm?Ecology 74(6):1659–1673.doi:10.2307/1939924

Liski J,Lehtonen A,Palosuo T et al(2006)Carbon accumulation in Finland’s forests 1922–2004-an estimate obtained by combination of forest inventory data with modelling of biomass,litter and soil.Ann For Sci 63(7):687–697.doi:10.1051/forest:2006049

Liu C,Zhang L,Li F et al(2014)Spatial modeling of the carbon stock of forest trees in Heilongjiang Province,China.J For Res 25(2):269–280.doi:10.1007/s11676-014-0458-x

Matheron G (1963) Principles of geostatistics.Econ Geol 58(8):1246–1266.doi:10.2113/gsecongeo.58.8.1246

McDonald RI,Peet RK,Urban DL(2003)Spatial pattern of Quercus regeneration limitation and Acer rubrum invasion in a Piedmont forest. J Veg Sci 14(3):441–450. doi:10.1658/1100-9233(2003)014%5B0441:SPOQRL%5D2.0.CO;2

Montes F,Herna′ndez MJ,Can?ellas I(2005)A geostatistical approach to cork production sampling estimation in Quercus suber forests.Can J For Res 35(12):2787–2796.doi:10.1139/x05-197

Mora A,Beer J(2013)Geostatistical modeling of the spatial variability of coffee fi ne roots under Erythrina shade trees and contrasting soil management.Agrofor Syst 87(2):365–376.doi:10.1007/s10457-012-9557-x

Moran PA (1950)Notes on continuous stochastic phenomena.Biometrika.doi:10.1093/biomet/37.1-2.17

Neilson ET,MacLean DA,Meng FR et al(2007)Spatial distribution of carbon in natural and managed stands in an industrial forest in New Brunswick,Canada.For Ecol Manag 253(1):148–160.doi:10.1016/j.foreco.2007.07.017

Nelson TA,Boots B(2008)Detecting spatial hot spots in landscape ecology.Ecography 31(5):556–566.doi:10.1111/j.0906-7590.2008.05548.x

Ord JK,Getis A(1995)Local spatial autocorrelation statistics:distributional issues and an application. Geogr Anal 27(4):286–306.doi:10.1111/j.1538-4632.1995.tb00912.x

Ord JK,Getis A(2001)Testing for local spatial autocorrelation in the presence of global autocorrelation.J Reg Sci 41(3):411–432.doi:10.1111/0022-4146.00224

Pan Y,Birdsey RA,Fang J et al(2011)A large and persistent carbon sink in the world’s forests.Science 333(6045):988–993.doi:10.1126/science.1201609

Rossi RE,Mulla DJ,Journel AG et al(1992)Geostatistical tools for modeling and interpreting ecological spatial dependence.Ecol Monogr 62(2):277–314.doi:10.2307/2937096

Saatchi SS,Harris NL,Brown S et al(2011)Benchmark map of forest carbon stocks in tropical regions across three continents.Proc Natl Acad Sci USA 108(24):9899–9904.doi:10.1073/pnas.1019576108

Sales MH,Souza CM,Kyriakidis PC et al(2007)Improving spatial distribution estimation of forest biomass with geostatistics:a case study for Rondo?nia,Brazil.Ecol Model 205(1):221–230.doi:10.1016/j.ecol model.2007.02.033

Sarangi A,Cox CA,Madramootoo CA(2005)Geostatistical methods for prediction of spatial variability of rainfall in a mountainous region.TASABE 48(3):943–954.doi:10.13031/2013.18507

Tang QY,Zhang CX(2013)Data processing system(DPS)software with experimental design,statistical analysis and data mining developed foruse in entomological research.InsectSci 20(2):254–260.doi:10.1111/j.1744-7917.2012.01519.x

Webster R,Oliver MA(2007)Geostatistics for environmental scientists.Wiley,New York

Zawadzki J,Cieszewski CJ,Zasada M et al(2005)Applying geostatistics for investigations of forest ecosystems using remote sensing imagery.Silva Fenn 39(4):599.doi:10.14214/sf.369

Zhang L,Gove JH,Heath LS(2005)Spatial residual analysis of six modeling techniques.Ecol Model 186(2):154–177.doi:10.1016/j.ecolmodel.2005.01.007

Zhang J,Ge Y,Chang J et al(2007)Carbon storage by ecological service forests in Zhejiang Province,subtropical China.For Ecol Manag 245(1):64–75.doi:10.1016/j.foreco.2007.03.042

Zhang L,Ma Z,Guo L(2009)An evaluation of spatial autocorrelation and heterogeneity in the residuals of six regression models.For Sci 55(6):533–548

Zhen Z,Li F,Liu Z et al(2013)Geographically local modeling of occurrence,count,and volume of downwood in Northeast China.Appl Geogr 37:114–126.doi:10.1016/j.apgeog.2012.11.003