• <tr id="yyy80"></tr>
  • <sup id="yyy80"></sup>
  • <tfoot id="yyy80"><noscript id="yyy80"></noscript></tfoot>
  • 99热精品在线国产_美女午夜性视频免费_国产精品国产高清国产av_av欧美777_自拍偷自拍亚洲精品老妇_亚洲熟女精品中文字幕_www日本黄色视频网_国产精品野战在线观看 ?

    Spatial variability of soil organic carbon in the forestlands of northeast China

    2014-09-06 11:13:31LingLiuHaiyanWangWeiDaiXiangdongLeiXiaojuanYangXuLi
    Journal of Forestry Research 2014年4期

    Ling Liu · Haiyan Wang · Wei Dai · Xiangdong Lei Xiaojuan Yang · Xu Li

    ORIGINAL PAPER

    Spatial variability of soil organic carbon in the forestlands of northeast China

    Ling Liu · Haiyan Wang · Wei Dai · Xiangdong Lei Xiaojuan Yang · Xu Li

    Received: 2013-10-17; Accepted: 2014-03-06

    ? Northeast Forestry University and Springer-Verlag Berlin Heidelberg 2014

    Soil organic carbon (SOC) is an effective indicator of soil fertility and productivity, and it varies spatially and temporally in relation to other soil properties. Spatial variability of SOC in the forestlands of northeast China was characterized using geostatistics. Soil samples at the depths of 0-20 cm, 20-40 cm and 40-60 cm were collected from sixty-three temporary plots to evaluate SOC concentration and density (SOCD) and other soil properties. We analyzed correlations between SOC and soil properties. Soil organic carbon concentrations were high. The total amount of C stored in soil (0-60 cm) was 16.23 kg·m-2with the highest SOCD of 7.98 kg·m-2in topsoil. Soil properties in most cases differed by horizon, suggesting different processes and effects in each horizon. Soil organic carbon had positive relationships with total N, P and K as well as readily available K, but did not show a significant positive correlation with available P. Spatial factors including elevation, slope and aspect affected SOC distribution. Soil organic carbon at 0-60 cm had strong spatial autocorrelation with nugget/sill ratio of 5.7%, and moderate structured dependence was found at 0-20 cm, which indicatedthe existence of a highly developed spatial structure. Spatial distributionsof SOC concentration and SOCD were estimated using regression-kriging, with higher prediction accuracy than ordinary kriging. The fractal dimension of SOC indicated the preferential pattern of SOC distribution, with the greatest spatial heterogeneity and strongest spatial dependence in the northeast-southwest direction.

    northeast China, soil organic carbon, spatial variability, geostatistics

    Introduction

    Soil organic carbon (SOC) is a key determinant of soil quality (Stocking 2003) and represents one of the largest reservoirs of organic carbon in the global carbon cycle that can influence global warming (Lützow et al. 2006). In terrestrial ecosystems, forest soil carbon is the largest carbon sink. Soil organic carbon storage in global forest ecosystems was estimated as 402×1012-787×1012kg, or 25%-50% of the global soil carbon storage (Arrouays and Pelissier 1994). In northeast China, considerable organic carbon is stored in soils because of the dense forest vegetation and the thick accumulation of organic matter. Humic cambisols (Dark brown forest soils, in Chinese soil taxonomy) are typical soils mainly found in warm-temperate, mixed coniferous and broad-leaved forests and cold-temperate coniferous forests in northeast China. The soils were higher in natural fertility and supported dense vegetation. The soils of 40.2×106ha are mainly located in the provinces of Heilongjiang, Jilin and Inner Mongolia (Zhao et al. 2004). These areas have been, therefore, among the most important for timber production and as a sink for organic matter on a global scale. Humic cambisols are also considered to store substantial amounts of organic matter, thereby serving as CO2sinks.

    Geostatistics is an effective method for study of spatial distribution characteristics and their variability. Geostatistical methods quantify spatial distributions based on the spatial scale of the study area, distance between detecting points, and spatial patterns of modeling semivariograms. They have been widely ap-plied to evaluate spatial correlation in soils and to analyze the spatial variability of soil characteristics, such as soil physical and chemical properties (Vieira et al. 2007; Zheng et al. 2009) and biological properties (Fromm et al. 1993). Similar to other soil properties, SOC concentration and stock are highly variable over space (horizontally and vertically) and time (Wigginton et al. 2000). Soil organic carbon density can vary due to intrinsic or extrinsic factors. Intrinsic variability is caused by natural variation (Cambardella et al. 1994) and extrinsic variability results from external factors imposed on sites (Rao and Wagenet 1985). Spatial correlation of SOC can be explored by geostatistical analysis. On a national scale, Lacelle established the soil landscape and SOC description database of Canada, and calculated that SOC storage at depths of 0-30 cm was 72.8×1012kg (Lacelle 1998). Bernoux estimated 0-30 cm SOC pool of 36.4(±3.4)×1012kg in Brazil with soil and vegetation maps (Bernoux et al. 2002). Soil carbon reserves in different soil horizons were studied in New Zealand (Scott et al. 2002) and Congo (Schwartz and Namri 2002). In China, with the aid of GIS spatial analysis, soil organic carbon density (SOCD) and storage were estimated for the entire country and a high degree of spatial variability was documented. Many studies investigated the spatial variability of SOC on a regional scale. Xue et al. (2003) surveyed the Qaidam basin and calculated the SOCD to be 6.11 kg·m-2. Using current 1:1000000 digital soil maps and 736 typical soil profiles in northeast China, organic carbon storage in topsoil and its spatial distribution were estimated, and the average SOCD at 0-1 m depths was 16.13 kg·m-2(Sun et al. 2004). In the tropical area of south China, SOC storage and spatial distribution were studied by combining soil data with a digital elevation model (DEM). Similar studies were conducted in Jiangsu, Guangdong and Anhui provinces (Jiang et al. 2005; Gan et al. 2003; Cheng and Xie 2009). According to previous studies, the initial application of GIS was to soil and land resource survey and evaluation, with focus on assessment of farmland soil degradation and mapping. Forestland, as a high potential land use type for carbon storage, has not yet gained much attention. Moreover, forest SOC storage or SOCD estimation on regional scales has large uncertainties because it is mainly based on soil survey data and literature and lacks detailed and reliable measurement. Ordinary kriging (OK) has been widely used for spatial interpolation in plains, but its predictive capability is limited in regions of complex topography with highly developed spatial conditions such as mountains and hills. Some ancillary variables such as terrain attributes and vegetation indices extracted from DEMs and remote sensing images are available for digital soil maps. Combining maps of soil types with soil properties information from point observations can improve spatial prediction of soil properties compared with prediction from a soil map only (Goovaerts and Journel 1995).

    The forest area in northeast China covers 3590×106ha, accounting for 24% of the total forest area of China, and the total forest stock volume amounts to 3.2×109m3(Department of Forest Resources Management 2010). It is the largest forest area and the transition region between broad-leaved mixed forest to the south and pure coniferous forest to the north. Research in this region has mainly focused on the relationship between forest and climate change and on spatial heterogeneity of organic carbon storage in vegetation layer. However, the forest soil carbon pool has seldom been studied and most spatial interpolation used OK. Information on the spatial distribution of SOC is important for evaluating soil ecological functions and understanding soil carbon sequestration processes. The use of geostatistics in combination with DEM and regression-kriging (RK) can provide such information and help to define and improve the accuracy of representations of the distribution of soil parameters across landscapes. Increased accuracy in the delimitation of SOC content classes could assist in meaningful soil quality prediction and assessment of the impacts of soil redistribution. The objectives of the present study were to evaluate the carbon-related soil properties and terrain attributes, and to determine the spatial variability and SOC distributions in different soil horizons of Humic cambisols in forests of northeast China using geostatistical analysis. This information will provide a scientific basis for sustainable forest management and ecological environment improvement in this region.

    Materials and methods

    Study area

    The study was carried out in Jincang Forest Farm, a typical forest area in the north temperate region of Wangqing County, Jilin Province, China (Fig. 1). The study area was in the eastern part of the low-middle hill region of the Changbai Mountain range (43°19′-43°23′ N, 130°26′-130°37′ E) at elevations of 562 to 950 m. The study area has a north temperate, monsoon climate with hot summers, cold winters, and abrupt changes of temperature. Average annual air temperature is 3.9°C. Average annual precipitation is 547 mm, concentrated from July to September (1970–2010). Soils are characterized by thick and black topsoil, loamy texture, an acidic pH range, large amount of humus and high natural fertility. The predominant soils in the area are Humic cambisols, together with meadow, peaty, swamp soil in some valleys that are developed from residual deposits of parent material such as granite, gneiss and basalt. Due to special landforms, abundant precipitation and changeable climate, rich vegetation resources are formed, including coniferous forest, broad-leaved forest, and mixed coniferous and broad-leaved forest. The major arbor species are Changbai larch (Larix olgensis), spruce (Picea koraiensis), Korean pine (Pinus koraiensis Sieb. et Zucc.), maple (Acer mono Maxim.), ash (Fraxinus mandshurica Rupr.), Asian white birch (Betula platyphylla Suk.) and linden (Tilia amurensis Rupr.).

    Soil sampling and analysis

    In August 2011, a total of 63 temporary circular sample plots (250 m2each plot) were set up in the study area. Within each plot, a soil profile was dug to describe profile characteristics and samples from the “bottom up” were taken with a ring sampler formeasurement of soil bulk density. The surface litter was removed prior to composite soil sampling. Soil samples at depths of 0-20, 20-40, and 40-60 cm were collected with 5 soil cores each well mixed into a composite soil sample, which was further divided into two sets of sub-samples. One set of the sub-samples was air dried and passed through a 2 mm sieve for analysis of soil texture, pH and available P and K; the remaining sub-sample was air-dried, homogenized and passed through a 0.25 mm sieve for analyses of SOC, total N, P and K (Li et al. 2011). Gravel and stones (fractions >2 mm) were sorted out and weighed. Soil texture was measured using the hydrometer method (Wilde et al. 1972). Soil pH was measured (1:2.5 = w/v) with a pH-meter (pH/ISE meter, model 710A; Orion). Soil organic carbon (SOC) concentration was determined by Walkley-Black wet oxidation method (Nelson and Sommers 1996; Bao 2000). Soil total N (TN) concentration was determined by Kjeldahl’s azotometer and total P (TP) by spectrophotometer, following wet digestion in concentrated H2SO4(Bremner 1996). Total K (TK) was determined by flame spectrophotometer following wet digestion in HF-HClO4(Jackson 1958; Knudsen et al. 1982). Available P (AP) was extracted with HCl-NH4F and estimated by spectrophotometer (Bray and Kurtz 1945). Readily available K (AK) was extracted with 1.0 M NH4OAc and quantified by flame spectrophotometer (Carson 1980). Soil bulk density was determined by oven-drying the soil samples collected with a ring sampler at 105°C for more than 8 h. In addition to soil sampling, the longitudes and latitudes of the sampling points were recorded with global positioning system (GPS). Forest type and survey, including environmental information, vegetation and management, for each sampling point were also recorded. A digital elevation model (DEM) of the study area was acquired at 5-m spatial resolution. The DEM and sampling points were geo-referenced to the Gauss Kruger coordinate system with Pulkovo 1942 GK Zone 22 (Fig. 1).

    Fig. 1 Study area location and sampling design.

    Calculation of soil organic carbon density (SOCD)

    Soil organic carbon density (SOCD) is the organic carbon content per unit area. We calculated SOCD in the profile with k layers (Sk, kg·m-2), on a ground-area basis up to 60 cm depth as follows:

    where, Ciis SOC concentration (g·kg-1), Diis soil bulk density (g·cm-3), Eiis soil thickness (cm) and Giis volumetric percentage of the fraction >2 mm (rock fragments). Siis soil organic carbon density in the ith horizon. In the present study, k=3.

    Statistical analysis

    Descriptive statistics

    A descriptive statistical analysis was carried out to determine the mean, minimum, maximum, standard deviation (SD) and coefficient of variation (CV) of the original data variables. The correlations among SOC, other soil properties and topographic factors such as elevation and aspect were also tested using correlation analysis and multiple linear stepwise regression. The significance of the difference between means of each property variable for different soil horizons was analyzed using one-way ANOVA. The statistical software SPSS 18.0 was used for statistical analysis.

    Geostatistical analysis

    The terrain analysis was based on a 5-m resolution digital elevation model (DEM) generated from a 1: 5000 contour topographical map. The elevation, slope in degrees, and aspect expressed in positive degrees from 0 to 360 measured clockwise from north of the study site were generated from the DEM using ArcGIS 9.3 (Esri Company, RedLands, USA). Geostatistical analysis was performed using the software GS+ and ArcGIS 9.3. A sample semivariogram was calculated for each soil variable as follows (Goovaerts 1999; Wang and Shao 2013):

    where, z(xi) and z(xi+h) are sample values at two points separated by the distance interval h, and N(h) is the total number of sample pairs for the lag interval h. For irregular sampling, it is rare for the distance between the sample pairs to be exactly equal to h. Therefore, h is often represented by a distance interval. γ(h) is the sample semivariance for the distance lag h, which reflects spatial relationship between neighbor sampling points.

    The sample semivariogram was calculated and then fitted with a suitable theoretical model chosen from spherical, exponential, Gaussian, Mate′rn or Stein’s Mate′rn (Cressie 1993). The semivariogram models provide information about the spatial structure as well as the input parameters for kriging interpolation. Spatial dependence of soil variables can be evaluated by three important parameters: nugget (C0), partial sill (C1) and range in the semivariogram plot. The nugget is defined as the variability at a scale smaller than the sampling interval and/or sampling and analyticalerror; the partial sill reflects the amount of spatial structural variance; and the range expresses the distance at which the semivariogram stabilizes around a limiting value. The spatial variation distance increasing to a constant value (sill) is known as range of spatial dependence (Cambardella et al. 1994), and it is affected by the sampling scale. The distribution distances of sampling points closer to the range are spatially correlated and the distances greater than the range are independent. In addition, the fractal dimension could be used to determine the spatial variability of soil properties. When the fractal dimension values are small, soil properties have weak spatial dependence and simple distribution, or vice versa (Burrough 1981, 1983; Wang and Shao 2013).

    Regression-kriging (RK) uses traditional regression methods that model the trend or drift in conjunction with geostatistical methods that model the random residuals from the regression (Hengl et al. 2004; Rivero et al. 2007). In RK the trend and residuals from the regression are fitted separately and then summed afterwards. It considers influence of geographic factors in the complex topography and simulates both spatial distribution trend and uncertainty. The trend item was separated by establishing regression equations between the auxiliary variable and target variable, and interpolation was applied to the residual with ordinary kriging (OK). The predicted value of the target variable was then obtained by spatial overlay between the trend item of regression prediction and OK value of the residual.

    Accuracy evaluation

    The prediction accuracy could be tested by comparing observed values with predicted values for soil nutrients, soil organic carbon (SOC) and soil organic carbon density (SOCD) from sample plots. Mean prediction error (M) and root mean square error (R) of prediction were adopted as the evaluation indices, and R was used to analyze the relative accuracy improvement value by comparing regression-kriging (RK) with ordinary kriging (OK) (Lal 2000; Kucharik et al. 2001).

    where, n is a sample point of test set, Zoiand Zpiare, respectively, observed values and predicted values, ROKis root mean square error of prediction using OK, and RRKis the prediction using RK. The prediction accuracy of RK is higher than OK when the R1value is positive, and the larger the value, the greater the accuracy. In this study, 10 soil sampling points were randomly extracted from 63 points to test the prediction accuracy of the remaining 53 points.

    Results

    Descriptive statistics of soil organic carbon (SOC) and other soil properties

    Table 1 shows the descriptive statistics of soil organic carbon (SOC) concentration and soil organic carbon density (SOCD). The mean of soil characteristics indicated the central tendency, and standard deviation (SD), coefficient of variation (CV) and skewness were used as the estimates of variability from each site. The average SOC concentration was 42.32 g·kg-1in the topsoil, and it decreased with depth until reaching 14.35 g·kg-1at 40-60 cm. The CV was low for SOC concentration in 0-20 cm (CV<10%), indicating low variability of SOC in topsoil. However, at soil depths of 20-40cm and 40-60 cm, SOC concentration showed high CVs (67% and 80%). Accordingly, the SOCD was 7.98 kg·m-2for the topsoil, and it declined gradually with increasing depth to 5.23 kg·m-2at 20-40 cm and 3.06 kg·m-2at 40-60 cm. Although there was a decreasing trend with increasing depth, the subsoil also contained a considerable amount of C, and the total amount of organic C stored to the depth of 60 cm was 16.23 kg·m-2, with a CV of 52%.

    Table 1: Descriptive statistics of soil organic carbon (n=63)

    The descriptive statistics for other soil properties in the study area are shown in Table 2. Soil pH had a mean of 5.41-5.57 with similar SDs of 0.26-0.30 in three soil layers. The CVs of total N (TN), total P (TP) and available P (AP) all exceeded 50%, which indicated considerable spatial variation in these properties within the study area. The concentrations of TN, TP and AP had means of 1.07 g·kg-1, 0.38 g·kg-1and 25.49 mg·kg-1, respectively, at 0-60 cm depth for all sampling points, and means and SDs varied by soil horizon. Mean values for all of these properties peaked at 0-20 cm and declined with increasing depth.

    Soil properties other than AP varied significantly by soil horizon. Concentrations of total N (TN), total P (TP), total K (TK) and readily available K (AK) declined with increasing depth in the following rank order: 0-20 cm > 20-40 cm > 40-60 cm withsignificant differences between each horizon (p<0.05); this trend was, however, not followed by available P (AP). In general, the typical distribution of soil chemical properties suggested higher variation between horizons and better nutrient status in the topsoil.

    Table 2: Descriptive statistics of soil properties in different soil horizons

    Correlation analysis The correlation matrix of selected soil properties was shown in Fig. 2. Soil organic carbon (SOC) concentration was significantly and positively correlated with total N (TN), total P (TP), total K (TK) and readily available K (AK), but not with available P (AP). Soil TN was positively related with SOC concentration (r = 0.795, p<0.01), but its correlation with AP was weak. Correlations of AP with TN and TP were weak, while correlation of AP with AK became stronger with increasing depth. In each horizon, there were strong positive correlations between soil properties. In addition, TN, TP, TK and AK had significant positive correlations with each other. The difference in SOC concentration was caused by natural organic matter distribution. To a certain extent, soil nitrogen levels affect SOC concentration through soil microorganisms. Different concentrations of P and K can be caused by different effects of organic matter, different pH values, and mineral composition.

    To the depth of 60 cm, the significant correlations between soil organic carbon (SOC) concentration and topographic attributes indicated that SOC concentration tended to be greater at higher elevation. The correlation between SOC concentration and elevation was very significant and positive (r =0.532, p<0.01). A very significant positive correlation was also found between soil organic carbon density (SOCD) and elevation (r = 0.565, p<0.01), but there was a very significant negative correlation between SOCD and slope (r = -0.451, p<0.01).

    Fig. 2: Correlation matrix between soil properties.??significant at probability level of 0.01;?significant at probability level of 0.05. TN is total N; TP is total P; TK is total K; AP is available P; AK is readily available K; SOC is soil organic carbon.

    There were very significant correlations between SOC concentration and SOCD (Table 3), and between these variables and cosine of slope aspect (p<0.01). SOC had, however, weak correlation with sine of slope aspect (which showed the degree of facing north), indicating that the more the aspect faced north the larger the SOC concentration. It tended to be greater on shady slopes at higher elevations. SOC in the topsoil had the same correlations with elevation and cosine of aspect, but this was not the case at depths of 20-40 cm and 40-60 cm.

    Table 3: Correlations between SOC concentration as well as SOCD and spatial factors (n=63)

    level of 0.05. SOC is soil organic carbon; SOCD is soil organic carbon density. cosα=cosine of aspect; sinα=sine of aspect. Aspect represents the degree of facing north.

    Multiple linear stepwise regression of soil organic carbon (SOC) concentration and soil organic carbon density (SOCD) on topographic factors yielded the following equations. Elevation and cosine of aspect (cosα) emerged as the optimal predictors of SOC and SOCD.

    Equations for 0-60 cm were:

    where, SOCis soil organic carbon concentration, SOCDis soil organic carbon density, and H is elevation. The coefficients of determination for fitting SOC and SOCD were 0.344 and 0.385, respectively, p<0.01 in both cases indicating good fits for both equations (6) and (7).

    Equations for 0-20 cm were:

    Geostatistical analysis

    Spatial variability of soil organic carbon (SOC) and semivariogram models

    Geostatistical parameters of soil organic carbon (SOC) are shown in Table 4. After interpolating the residual with ordinary kriging (OK), the log-transformed data of SOC concentration at 0-20 cm depth was fitted well by a Gaussian model having a nugget/sill ratio (NSR) of 5.7%, a range of 0.294 km, and very strong spatial autocorrelation. The large spatial range indicated a strongly structured regional pattern of SOC. The Gaussian model was best fitted to sample semivariograms of SOC at 0-60 cm depth with moderate structured dependence (NSR=47.8%). As for soil organic carbon density (SOCD), the NSRat 0-20 cm was 15.4%, suggesting a considerable development of the spatial structure and high spatial correlativity. SOC at 0-60 cm depth was fitted by an exponential model and the spatial structure was moderate to well developed (NSR=49.3%).

    Table 4: Semivariogram models for SOC concentration and SOCD’s residual error

    Accuracy evaluation for different prediction methods

    The regression-kriging (RK) method considers environmental factors and separates residual error to increase prediction accuracy. Fifty-three sampling points were randomly extracted from 63 points to conduct ordinary kriging (OK) interpolation for regression prediction residual error of soil organic carbon (SOC). Fig. 3 shows the spatial distribution of SOC concentration and soil organic carbon density (SOCD) interpolated using RK. At 0-60 cm, the spatial distribution of SOC mainly concentrated on northerly aspects. Taking the digital elevation model (DEM) in comparison, SOC concentration and SOCD were relatively rich at high elevations. Soil organic carbon at 0-20 cm showed that spatial distribution of SOC concentration had a number of high-value centers but was mainly concentrated in high elevation regions. The results suggest that slopes of northerly aspect served as sinks of organic matter and became the largest sources of SOC.

    Fig. 3: Spatial distributions of soil organic carbon (SOC) and soil organic carbon density (SOCD) interpolated by regression-kriging.

    The remaining 10 sampling points were used to test the accuracy of the two prediction methods. Regression-kriging (RK) prediction was much more detailed on spatial variation and more closely approximated spatial distribution condition of the actual value (Table 5). Positive R1values indicated higher accuracy using the RK method.

    Table 5: Comparison of accuracy using different prediction methods

    The fractal dimension values of soil organic carbon (SOC) in different soil horizons are shown in Table 6. The greatest isotropic fractal dimension was found at 0-20 cm, indicating the most complex spatial pattern of SOC and the highest degree of uniformity. The fractal dimension at 20-40 cm was the smallest, indicating a relatively simple spatial pattern and stronger spatial dependence than other soil horizons.

    In terms of slope aspect, the greatest fractal dimension of soil organic carbon (SOC) at 0-20 cm was found in the northwest-southeast (NW-SE) direction, and the smallest in the northeast-southwest (NE-SW) direction, which had the greatest spatial heterogeneity and the strongest spatial dependence, and was also the direction of the longest transect through the study site (Fig. 3). As for the coefficient of determination, NE-SW direction showed a preferential pattern of SOC distribution. At 20-40 cm, the NE-SW direction had the smallest fractal dimension with the strongest spatial dependence; east-west (E-W) direction was the preferential pattern distribution for SOC. The fractal dimension of SOC at 40-60 cm depth was greatest in the NW-SE direction; and smallest in the NE-SW direction, which had greatest spatial heterogeneity. The coefficient of determination showed that NE-SW had the preferential pattern of SOC distribution.

    Table 6: Fractal dimension of SOC in different soil layers

    Discussion

    The average soil organic carbon (SOC) concentration was 27.40 g·kg-1and soil organic carbon density (SOCD) to the depth of 60 cm was 16.23 kg·m-2, which is in agreement with previous studies that estimated SOCD in a range of 16.55-27.02 kg·m-2in this region (Zu et al. 2011; Sun et al. 2004). Because of differences inthe locations of study areas and management activities, SOC varies spatially. In addition, research methods are not standardized, and estimates are not free of uncertainty due to sampling depth and quantity.

    Soil organic carbon (SOC) and soil organic carbon density (SOCD) were higher in topsoil and declined with increasing depth. Distribution of SOC is related to vegetation types and succession, and human disturbance, so organic carbon in topsoil is greatly influenced by these factors. In the forest ecosystem, soil organic matter is mainly from litter and root decomposition. That surface soil can obtain more forest litter leads to the enrichment of SOC. The SOCD is associated with soil bulk density and stoniness (diameters >2 mm). Several studies found that SOCD increased with increasing soil depth (Yang and Wang 2005). However, some studies in forests suggest that the changes in SOC stocks are due to the development of vegetation restoration and succession. For instance, vegetation restoration can increase non-active organic carbon concentration in surface soil (Su et al. 2005). Succession of restored vegetation can also increase SOC and nitrogen concentration (Christian et al. 2001).

    Correlation analysis quantifies the tendency of variables such as SOC, soil properties, and terrain attributes to vary in similar ways. Soil properties have indirect influence on SOC by affecting litter decomposition, microbial activity and organic carbon mineralization rates (Fayez and Esmael 2006). Soil total N (TN) concentration can promote the accumulation of biomass in understory vegetation and litter generation, thereby facilitating SOC accumulation. As a result, SOC, elevation, aspect, TN and total K (TK) were all interrelated, suggesting that the geographical characteristics and soil properties were correlated with SOC. Soil organic carbon (SOC) in all soil horizons had a positive relationship with TN, total P (TP), TK and readily available K (AK), in agreement with previous observations (Mládková et al. 2005). A significant positive correlation between SOC and elevation indicated that SOC concentration tended to be greater at higher elevation. This is inconsistent with studies reporting greater SOC concentrations on the valley floor because of the deposition of water, solutes and sediments down along the slope (Mueller and Pierce 2003; Sumfleth and Duttmann 2008; Wang et al. 2001). Our findings agree with those of Tsuia et al. (2004) and Zhang et al. (2012) who documented increasing SOC with increasing elevation.

    Spatial patterns of soil organic carbon (SOC) suggest that organic matter dynamics in the field were affected by topography. The spatial variability of SOC and soil organic carbon density (SOCD) differed by soil horizon. The spatial distribution of surface SOC demonstrated more developed spatial structure than the other two horizons, implying that topsoil might be more sensitive to environmental disturbances, whereas spatial variability of SOC at depths of 20-40 cm and 40-60 cm was greatly affected by random factors and suggested considerably weaker spatial correlativity. By combining the kriging map with the digital elevation model (DEM), the SOC concentration was related to elevation, indicating greater SOC in the higher terrain area. Similar findings with regard to the importance of elevation to SOC were reported by Geng et al (2011). However, some research yielded the opposite results (Cheng et al. 2009; Zhang et al. 2012). This study showed that predictive capability of ordinary kriging (OK) was sometimes limited in regions of complex topography affected by intense cultivation and land use change. Recently, more studies have demonstrated regression-kriging (RK) can provide secondary information to improve the prediction accuracy (Lark and Webster 2006; Liu et al. 2010). It is possible that other useful ancillary information, such as soil and forest types, can be added in using RK to produce more accurate representations of the spatial distribution of SOC (Hengl et al. 2004). The RK method with elevation as the predictor gave more details in space than the OK method in predicting spatial distributions of SOC variables. In the present study, it was difficult to use the spatial information reliably to establish correlations with the SOC variables by using RK analysis because study area spanned a small range of elevations.

    Conclusions

    Soil organic carbon (SOC) concentration and soil organic carbon density (SOCD) decreased significantly with increasing soil depth. Mean SOC concentration was 27.40 g·kg-1and SOCD in the whole soil profile (0-60 cm) was 16.23 kg·m-2. Correlations between SOC and soil chemical properties varied by soil depth. Soil organic carbon concentration showed significant positive relationships with total N (TN), total P (TP), total K (TK) and readily available K (AK), but not with available P (AP). A significant positive correlation was found between SOC and elevation, as well as between SOC and slope aspect. Soil organic carbon tended to be higher on sunny slopes. In the surface horizon, spatial variation produced by random factors was larger and had very small spatial correlation. With increasing depth, SOC was less sensitive to external influence with smaller spatial variation. The effect of terrain factors on SOC was very complex. Resulting maps showed the spatial distribution of SOC in the study area based on the data of spatial dependency, and they can be used to delineate the areas with different impact of elevation.

    Arrouays H, Pelissier P. 1994. Modeling carbon storage profiles in temperate forest humic loamy soils of France. Soil Science, 157(3): 185-192.

    Bao SD. 2000. Soil and agricultural chemistry analysis. Beijing: China Agriculture Press, p. 495. (In Chinese)

    Bernoux M, Conceicao SCM, Volkoff B. 2002. Brazil’s soil carbon stock. Soil Science, 66: 888-896.

    Bray RH, Kurtz LT. 1945. Determination of total, organic and available forms of P in soils. Soil Science, 59: 39-45.

    Bremner JM. 1996. Nitrogen—total. In: Sparks DL (ed), Methods of Soil Analysis. Madison, WI: SSSA-ASA, pp. 1085-1121.

    Burrough P. 1981. Fractal dimensions of landscapes and other environmental data. Nature, 294: 240-242.

    Burrough P. 1983. Multiscale sources of spatial variation in soil. I. The application of fractal concepts to nested levels of soil variation. European Journal of Soil Science, 34: 577-597.

    Cambardella CA, Moorman TB, Novak JM, Parkin TB, Karlen DL, Turco RF, Konopka AE. 1994. Field-scale variability of soil properties in central Iowa soils. Soil Science Society of America Journal, 58: 1501-1511.

    Carson PL. 1980. Recommended potassium test. In: Dahnke WC (ed) Recommended chemical soil test procedures for the North Central region. North Dakota Agricultural Experimental Station Bulletin 499: 17-18.

    Cheng XF, Xie Y. 2009. Spatial distribution of soil organic carbon density in Anhui Province based on GIS. Sientia Geographica Sinica, 29(4): 540-544. (In Chinese)

    Christian PG, Michael GR, Robert MH. 2001. Tree species and soil textural controls on carbon and nitrogen mineralization rates. Soil Science Society of America Journal, 65: 1272-1279.

    Cressie NAC. 1993. Statistics for Spatial Data. New York: John Wiley & Sons p. 900.

    Department of Forest Resources Management. 2010. The 7th national forest resources inventory and the status of forest resources. Forest Resource Management, 2: 1-8. (In Chinese)

    Fayez R, Esmael A. 2006. Soil microbial activity and litter turnover in native grazed and ungrazed rangelands in a semiarid ecosystem. Biology and Fertility of Soils, 43: 76-82.

    Fromm H, Winter K, Filer J. 1993. The influence of soil type and cultivation system on the spatial distributions of the soil fauna and microorganisms and their interactions. Geoderma, 60: 109-118.

    Gan HH, Wu SH, Fan XD. 2003. Reserves and spatial distribution characteristics of soil organic carbon in Guangdong Province. Chinese Journal of Applied Ecology, 14(9): 1499-1502. (In Chinese)

    Geng GP, Gao P, Lü SQ, Zhang J. 2011. Spatial distribution of soil organic matter and total nitrogen in Matiyu small watershed in hilly area of middle southern Shandong Province. Science of Soil and Water Conservation, 9(6): 99-105. (In Chinese)

    Goovaerts P, Journel AG. 1995. Integrating soil map information in modeling the spatial variation of continuous soil properties. Soil Science, 46(3): 396-414.

    Goovaerts P. 1999 Geostatistics in soil science: state-of-the-art and perspectives. Geoderma, 89: 1-45.

    Hengl T, Heuvelink GBM, Stein A. 2004. A generic framework for spatial prediction of soil variables based on regression-kriging. Geoderma, 120: 75-93.

    Jackson ML. 1958. Soil chemical analysis. NJ, USA:Prentice Hall, Inc. Englewood Cliffs, p. 498.

    Jackson RB, Caldwell MM. 1993. The scale of nutrient heterogeneity around individual plants and its quantification with geostatistics. Ecology, 74(2): 612-614.

    Jiang XS, Pan JJ, Li XL. 2005. Organic carbon density and storage in topsoil of Jiangsu and spatial analysis. Chinese Journal of Soil Science, 36(8): 501-503. (In Chinese)

    Knudsen D, Peterson GA, Pratt P. 1982. Lithium, sodium, and potassium. In: Page AL (eds), Methods of Soil Analysis, Vol. 2. Madison, Wis:American Society of Agronomy, Soil Science Society of America, pp. 225-246.

    Kucharik CJ, Brye KR, Norman JM, Foley JA, Gower ST, Bundy LG. 2001. Measurements and modeling of carbon and nitrogen cycling in agroecosystems of southern Wisconsin: potential for SOC sequestration during the next 50 years. Ecosystems, 4: 237-258.

    Lacelle B. 1998. Canada’s soil organic carbon database. In: Lal R, Kimble JM, Follett RF, Stewart BA (eds), Advances in Soil Science: soil processes and the carbon cycle. CRC Press, pp. 93-101.

    Lal R. 2000. Carbon sequestration in dry lands. Annual of Arid Zone, 39(1): 1-10.

    Lark RM, Webster R. 2006. Geostatistical mapping of geomorphic variables in the presence of trend. Earth Surf Process Landforms, 31: 862-874.

    Li LJ, Zeng DH, Yu ZY, Fan ZP, Yang D, Liu YX. 2011. Impact of litter quality and soil nutrient availability on leaf decomposition rate in a semi-arid grassland of Northeast China. Journal of Arid Environments, 75:787-792.

    Liu DW, Wang ZM, Zhang B. 2006. Spatial distribution of soil organic carbon and analysis of related factors in croplands of the black soil region, Northeast China. Agricultural Ecosystem Environment, 113: 73-81. (In Chinese)

    Liu SL, Li Y, Wu JS. 2010. Spatial variability of soil microbial biomass carbon, nitrogen and phosphorus in a hilly red soil landscape in subtropical China. Soil Science and Plant Nutrition, 56:693-704. (In Chinese)

    Lützow Mv, K?gel-Knabner I, Ekschmitt K, Matzner E, Guggenberger G, Marschner B, Flessa H. 2006. Stabilization of organic matter in temperate soils: mechanisms and their relevance under different soil conditions-a review. European Journal of Soil Science, 57(4):426-445.

    Mládková L, Bor?vka L, Drábek O. 2005. Soil properties and toxic aluminium forms in acid forest soils as influenced by the type of vegetation cover. Soil Science and Plant Nutrition, 51:741-744.

    Mueller TG, Pierce FJ. 2003. Soil carbon maps: enhancing spatial estimates with simple terrain attributes at multiple scales. Soil Science Society of America Journal, 67: 258-267.

    Nelson DW, Sommers LE. 1996. Total carbon, organic carbon, and organic matter. In: Sparks DL, Page AL, etc. (eds), Methods of Soil Analysis, Part 3. Chemical methods. Wisconsin, WI, USA: Soil Science Society of America Book Series, Vol. 5. pp. 961-1010.

    Rao PSC, Wagenet RJ. 1985. Spatial variability of pesticides in field soils: methods for data analysis and consequences. Weed Science, 33(suppl. 2):18-24.

    Rivero RG, Grunwald S, Bruland GL. 2007. Incorporation of spectral data into multivariate geostatistical models to map soil phosphorous variability in a Florida wetland. Geoderma, 140: 428-433.

    Schwartz D, Namri M. 2002. Mapping the total organic carbon in the soils of the Congo. Global and Planetary Change, 33: 77-93.

    Scott NA, Tare KR, Giltrap DJ. 2002. Monitoring land-use effects on soil carbon in New Zealand: quantifying baseline soil carbon stocks. Environmental Pollution, 116: 167-186.

    Stocking MA. 2003. Tropical soils and food security: the next 50 years. Science, 302(5649): 1356-1359.

    Su J, Zhao SW, Ma JD. 2005. Influence of man-made vegetation on carbon pool in southern Ningxia Region in Loess Plateau. Research of Soil and Water Conservation, 12(3): 50-52. (In Chinese)

    Sumfleth K, Duttmann R. 2008. Prediction of soil property distribution in paddy soil landscapes using terrain data and satellite information as indicators. Ecological Indicators, 8: 485-501.

    Sun WX, Shi XZ, Yu DS. 2004. Estimation of soil organic carbon density and storage of northeast China. Acta Pedologica Sinica, 41(2):298-300. (In Chinese)

    Tsuia C, Chena Z, Hsieh C. 2004. Relationships between soil properties and slope position in a lowland rain forest of southern Taiwan. Gerderma, 123: 131-142.

    Vieira VA, Mello CR, Lima JM. 2007. Spatial variability of soil physical attributes in small watershed. Ciência e Agrotecnologia, 31(5): 1477-1485. (In Portuguese).

    Wang J, Fu BJ, Qiu Y. 2001. Soil nutrients in relation to land use and landscape position in the semi-arid small catchment on the loess plateau in China. Journal of Arid Environments, 48(4): 537-550

    Wang YQ, Shao MA. 2013. Spatial variability of soil physical properties in a region of the Loess Plateau of PR China subject to wind and water erosion. Land Degradation and Development, 24(3): 296-304.

    Wigginton JD, Lockaby BG, Trettin CC, Nelson EA, Kolka RK, Wisniewski J. 2000. Soil organic matter formation and sequestration across a forested floodplain chronosequence. In: Proceedings of a workshop held at Clemson University, 12-14 April 1999, pp. 141-155.

    Wilde SA, Voigt GK, Iyer JC. 1972. Soil and plant analysis for tree culture. New Delhi, India: Oxford and IBH Publishing Co., p. 209.

    Xue L, Ma HZ, Cao GC, Sha ZJ, Ou LY. 2003. The application of geographic information systems on the calculation of storage of organic carbon in soils. Ecology & Environment, 12(4):419-422. (In Chinese)

    Yang JY, Wang CK. 2005. Soil carbon storage and flux of temperate forest ecosystems in northeastern China. Acta Ecologica Sinica, 25(11): 2875-2882. (In Chinese)

    Zhang W, Wang KL, Chen HS, He XY, Zhang JG. 2012. Ancillary information improves kriging on soil organic carbon data for a typical karst peak cluster depression landscape. Journal of the Science of Food and Agriculture, 92(5): 1094-1102.

    Zhao J, Zhang JM, Meng K. 2004. Spatial heterogeneity of soil nutrients in black soil, China-a case study at Hailun County. Bulletin of Soil and Water Conservation, 24(6): 53-57. (In Chinese)

    Zhao X, Wang Q, Kakubari Y. 2009. Stand-scale spatial patterns of soil microbial biomass in natural cold-temperate beech forests along an elevation gradient. Soil Biology and Biochemistry, 41: 1466-1474.

    Zheng Z, Zhang FR., Ma FY. 2009. Spatiotemporal changes in soil salinity in a drip irrigated field. Geoderma, 149: 243-248.

    Zu YG, Li R, Wang WJ. 2011. Soil organic and inorganic carbon contents in relation to soil physicochemical properties in northeastern China. Acta Ecologica Sinica, 31(18): 5207-5216. (In Chinese)

    DOI 10.1007/s11676-014-0533-3

    Project funding: The study was jointly supported by Natural Science Foundation of China (No. 31270697), the Fundamental Research Funds for the Central Universities (TD2011-2), State Forestry Administrative public service sector project “Key management techniques for the health of typical forest types in China” (20100400201) and National ‘973’project “Soil carbon stock and its temporal and spatial distribution pattern in natural forests” (2011CB403201).

    The online version is available at http:// www.springerlink.com

    Ling Liu, Haiyan Wang (), Wei Dai, Xiaojuan Yang, Xu Li

    College of Forestry, Beijing Forestry University, Beijing 100083, P. R. China. E-mail: haiyanwang72@aliyun.com

    Xiangdong Lei

    Institute of Forest Resource Information Techniques, Chinese Academy of Forestry, Beijing 100091, P. R. China.

    Corresponding editor: Zhu Hong

    啦啦啦观看免费观看视频高清| 国产黄片视频在线免费观看| 日本av手机在线免费观看| 日日干狠狠操夜夜爽| 亚洲国产精品成人综合色| 色播亚洲综合网| 国内精品宾馆在线| a级毛片免费高清观看在线播放| 久久中文看片网| av.在线天堂| 国产伦精品一区二区三区四那| 国产精品久久电影中文字幕| 亚洲七黄色美女视频| 亚洲国产欧美在线一区| 日韩精品有码人妻一区| 欧美日韩精品成人综合77777| 级片在线观看| 成人无遮挡网站| 国内揄拍国产精品人妻在线| 亚洲国产精品sss在线观看| 联通29元200g的流量卡| 直男gayav资源| 老司机福利观看| 亚洲成人久久爱视频| 五月伊人婷婷丁香| 日韩欧美精品v在线| 91狼人影院| 国产亚洲av片在线观看秒播厂 | 男人舔奶头视频| 大又大粗又爽又黄少妇毛片口| 国产人妻一区二区三区在| 亚洲欧洲日产国产| 热99在线观看视频| 国产亚洲精品久久久久久毛片| 国产精品久久电影中文字幕| 天美传媒精品一区二区| 亚州av有码| 色噜噜av男人的天堂激情| 欧美+亚洲+日韩+国产| 身体一侧抽搐| 一个人看视频在线观看www免费| 国产精品久久久久久久久免| 久久久色成人| 久久久精品94久久精品| 国产 一区精品| 欧美区成人在线视频| 成人鲁丝片一二三区免费| 亚洲精品久久国产高清桃花| 美女内射精品一级片tv| 69av精品久久久久久| 久久久久久九九精品二区国产| 少妇熟女aⅴ在线视频| 欧美日本视频| 久久人妻av系列| 国产真实乱freesex| 卡戴珊不雅视频在线播放| 亚洲欧美日韩卡通动漫| 长腿黑丝高跟| 国产亚洲精品久久久久久毛片| 成人特级黄色片久久久久久久| 亚洲人成网站在线观看播放| 老司机福利观看| 人妻少妇偷人精品九色| 成人无遮挡网站| 女人十人毛片免费观看3o分钟| 亚洲精品色激情综合| 久久精品国产清高在天天线| 精品人妻视频免费看| 国产爱豆传媒在线观看| 插逼视频在线观看| 国产精品美女特级片免费视频播放器| 99在线视频只有这里精品首页| 亚洲18禁久久av| 久久99精品国语久久久| 人妻制服诱惑在线中文字幕| 精品久久久噜噜| 综合色av麻豆| av天堂在线播放| 大香蕉久久网| 欧美激情国产日韩精品一区| 午夜福利在线观看吧| 色综合站精品国产| 欧美性感艳星| 久久热精品热| 夜夜看夜夜爽夜夜摸| 内地一区二区视频在线| 在线免费观看不下载黄p国产| 观看美女的网站| 黄色日韩在线| 看片在线看免费视频| 可以在线观看的亚洲视频| 少妇的逼好多水| 精华霜和精华液先用哪个| 熟女电影av网| 欧美3d第一页| 麻豆精品久久久久久蜜桃| 欧美成人一区二区免费高清观看| 最好的美女福利视频网| 午夜福利在线在线| 国产毛片a区久久久久| 久久久精品大字幕| 一本久久精品| 日本黄色视频三级网站网址| kizo精华| 高清在线视频一区二区三区 | 国产亚洲精品av在线| 亚洲av免费高清在线观看| 日本五十路高清| 边亲边吃奶的免费视频| 国产伦一二天堂av在线观看| 欧美激情久久久久久爽电影| 一级毛片aaaaaa免费看小| 国产91av在线免费观看| 免费搜索国产男女视频| 青春草亚洲视频在线观看| 97超碰精品成人国产| 啦啦啦观看免费观看视频高清| 国产精华一区二区三区| 色噜噜av男人的天堂激情| 午夜福利高清视频| 欧美极品一区二区三区四区| 人人妻人人澡人人爽人人夜夜 | 波多野结衣巨乳人妻| 日韩 亚洲 欧美在线| videossex国产| 国产精品精品国产色婷婷| 在线国产一区二区在线| 91精品国产九色| 嘟嘟电影网在线观看| 爱豆传媒免费全集在线观看| 亚洲国产欧洲综合997久久,| 国产私拍福利视频在线观看| 国产不卡一卡二| 伦理电影大哥的女人| 秋霞在线观看毛片| 春色校园在线视频观看| 一区二区三区四区激情视频 | 欧美潮喷喷水| 亚洲国产精品久久男人天堂| 亚洲av熟女| 亚洲av第一区精品v没综合| 午夜a级毛片| 精品不卡国产一区二区三区| 国内精品久久久久精免费| 久久人妻av系列| 日本三级黄在线观看| or卡值多少钱| 国国产精品蜜臀av免费| 亚洲精品亚洲一区二区| 国产在线精品亚洲第一网站| 又黄又爽又刺激的免费视频.| 日韩av在线大香蕉| 久久精品91蜜桃| 日韩一区二区视频免费看| 欧美一级a爱片免费观看看| 99精品在免费线老司机午夜| 久久国内精品自在自线图片| 免费看美女性在线毛片视频| 亚洲精品粉嫩美女一区| 亚洲一级一片aⅴ在线观看| 国产精品精品国产色婷婷| 丰满的人妻完整版| 波多野结衣高清作品| 久久亚洲精品不卡| 天天躁夜夜躁狠狠久久av| 久久久久久久久大av| 亚洲精品乱码久久久久久按摩| 午夜福利高清视频| 不卡一级毛片| 亚洲真实伦在线观看| 国产精品无大码| 久久精品国产自在天天线| 日韩三级伦理在线观看| av国产免费在线观看| 国产午夜精品久久久久久一区二区三区| 国产真实伦视频高清在线观看| 精品一区二区免费观看| 免费观看在线日韩| 2022亚洲国产成人精品| 小说图片视频综合网站| 国产爱豆传媒在线观看| 九九在线视频观看精品| 国产成人一区二区在线| 国产爱豆传媒在线观看| 国产三级中文精品| 亚洲经典国产精华液单| 欧美日韩一区二区视频在线观看视频在线 | av又黄又爽大尺度在线免费看 | 97在线视频观看| 免费av毛片视频| 美女国产视频在线观看| 国产成人福利小说| 日韩欧美 国产精品| 国产精品1区2区在线观看.| h日本视频在线播放| 欧美区成人在线视频| 国产精品99久久久久久久久| 欧美日韩一区二区视频在线观看视频在线 | 成人二区视频| 一级黄片播放器| 国产一级毛片在线| 又爽又黄无遮挡网站| 成人综合一区亚洲| 国产亚洲av嫩草精品影院| 免费观看人在逋| 一区二区三区免费毛片| 最近中文字幕高清免费大全6| 欧美激情在线99| 午夜亚洲福利在线播放| 亚洲精品影视一区二区三区av| 狂野欧美激情性xxxx在线观看| 精品午夜福利在线看| 日本爱情动作片www.在线观看| 99在线人妻在线中文字幕| 日本免费一区二区三区高清不卡| 成人特级黄色片久久久久久久| 日本熟妇午夜| 国产高清三级在线| 51国产日韩欧美| 99热这里只有精品一区| 国产v大片淫在线免费观看| 国产色婷婷99| 亚洲无线观看免费| 亚洲精品456在线播放app| 成年免费大片在线观看| 久久精品夜色国产| 久久久国产成人免费| а√天堂www在线а√下载| 免费观看人在逋| av在线播放精品| 国产精品免费一区二区三区在线| 欧美性感艳星| 免费搜索国产男女视频| av卡一久久| 日韩av在线大香蕉| 色播亚洲综合网| 日本撒尿小便嘘嘘汇集6| 有码 亚洲区| 亚洲精品色激情综合| 国产一区二区亚洲精品在线观看| 欧美zozozo另类| 天天躁日日操中文字幕| 亚洲高清免费不卡视频| av视频在线观看入口| 成人性生交大片免费视频hd| 免费搜索国产男女视频| 99在线视频只有这里精品首页| 欧洲精品卡2卡3卡4卡5卡区| 久久久欧美国产精品| 观看免费一级毛片| av黄色大香蕉| 免费av观看视频| 网址你懂的国产日韩在线| av在线亚洲专区| 你懂的网址亚洲精品在线观看 | 午夜精品在线福利| 国产v大片淫在线免费观看| 别揉我奶头 嗯啊视频| 午夜久久久久精精品| 成人美女网站在线观看视频| 成人午夜精彩视频在线观看| 欧美激情国产日韩精品一区| 看十八女毛片水多多多| 国产精品一区二区三区四区免费观看| 精品熟女少妇av免费看| 中文资源天堂在线| 国产伦精品一区二区三区四那| 国产精品电影一区二区三区| 日本在线视频免费播放| 男女那种视频在线观看| 精品少妇黑人巨大在线播放 | 麻豆乱淫一区二区| 99热这里只有精品一区| av在线亚洲专区| 国产午夜精品论理片| 小说图片视频综合网站| 亚洲精品乱码久久久久久按摩| 成人毛片a级毛片在线播放| 18禁黄网站禁片免费观看直播| 菩萨蛮人人尽说江南好唐韦庄 | 男女那种视频在线观看| 久久久久久久久久久丰满| 日韩欧美在线乱码| av专区在线播放| 日韩,欧美,国产一区二区三区 | 十八禁国产超污无遮挡网站| 国产精品电影一区二区三区| 天堂网av新在线| 波野结衣二区三区在线| 99国产极品粉嫩在线观看| 日韩欧美一区二区三区在线观看| 色吧在线观看| 欧美激情在线99| 国产人妻一区二区三区在| 99久久精品国产国产毛片| 99在线视频只有这里精品首页| 综合色丁香网| 日韩一区二区三区影片| 国产亚洲av片在线观看秒播厂 | 色视频www国产| 看片在线看免费视频| 两性午夜刺激爽爽歪歪视频在线观看| 美女脱内裤让男人舔精品视频 | 国国产精品蜜臀av免费| 久久久色成人| 国产单亲对白刺激| 成人午夜精彩视频在线观看| 国产精品99久久久久久久久| 在线观看av片永久免费下载| 亚洲欧洲国产日韩| 99久久精品国产国产毛片| 亚洲精品456在线播放app| 在线天堂最新版资源| 深夜精品福利| 国产中年淑女户外野战色| 国产精品不卡视频一区二区| 午夜福利高清视频| 亚洲无线在线观看| 国产黄色小视频在线观看| 日本在线视频免费播放| 欧美高清性xxxxhd video| 一本精品99久久精品77| 91狼人影院| 久久久久久久久久久免费av| av天堂在线播放| 亚州av有码| 国产一区二区激情短视频| 有码 亚洲区| 久久午夜福利片| 精品人妻一区二区三区麻豆| 欧美日韩综合久久久久久| 亚洲成人久久性| avwww免费| 欧美又色又爽又黄视频| 色哟哟哟哟哟哟| 少妇人妻一区二区三区视频| 日本免费a在线| 精品午夜福利在线看| 亚洲国产欧洲综合997久久,| 午夜免费激情av| 久久人人爽人人片av| 成人综合一区亚洲| 国产女主播在线喷水免费视频网站 | 午夜爱爱视频在线播放| 久久人人爽人人片av| 嘟嘟电影网在线观看| 麻豆成人av视频| 两性午夜刺激爽爽歪歪视频在线观看| 国产白丝娇喘喷水9色精品| ponron亚洲| 欧美日韩在线观看h| 国产三级中文精品| 国产男人的电影天堂91| 嫩草影院新地址| 国产亚洲精品久久久com| 亚洲国产精品久久男人天堂| 国产免费一级a男人的天堂| 黄片无遮挡物在线观看| АⅤ资源中文在线天堂| 亚洲欧美中文字幕日韩二区| 高清在线视频一区二区三区 | 亚洲国产高清在线一区二区三| 黄片wwwwww| 欧美又色又爽又黄视频| 一边亲一边摸免费视频| 综合色av麻豆| 在线免费十八禁| 国产av麻豆久久久久久久| 大香蕉久久网| 人人妻人人澡欧美一区二区| 亚洲一级一片aⅴ在线观看| 亚洲五月天丁香| 国产91av在线免费观看| 伦理电影大哥的女人| 国产精品久久久久久精品电影| 亚洲熟妇中文字幕五十中出| av卡一久久| 国产精品永久免费网站| 久久精品国产亚洲av涩爱 | 美女国产视频在线观看| 久久九九热精品免费| 亚洲第一电影网av| 亚洲中文字幕日韩| 只有这里有精品99| 久久久久久伊人网av| 中国美女看黄片| 国产高清三级在线| 直男gayav资源| 网址你懂的国产日韩在线| www.色视频.com| 嫩草影院精品99| 1024手机看黄色片| 热99re8久久精品国产| 亚洲人成网站在线观看播放| 日韩人妻高清精品专区| 中文精品一卡2卡3卡4更新| 69av精品久久久久久| 在线播放国产精品三级| 国产精品不卡视频一区二区| 亚洲欧洲国产日韩| 美女xxoo啪啪120秒动态图| 色视频www国产| 亚洲在线观看片| 在线a可以看的网站| 长腿黑丝高跟| 午夜爱爱视频在线播放| 亚洲人成网站在线播| 亚洲精华国产精华液的使用体验 | 天堂中文最新版在线下载 | 波野结衣二区三区在线| 国产精品一及| 午夜亚洲福利在线播放| 亚洲四区av| 国产片特级美女逼逼视频| 蜜桃亚洲精品一区二区三区| 51国产日韩欧美| 日韩视频在线欧美| 亚洲激情五月婷婷啪啪| 久久综合国产亚洲精品| 国产av一区在线观看免费| 国产成人午夜福利电影在线观看| 亚洲成av人片在线播放无| 国产精品无大码| 我要看日韩黄色一级片| 成人性生交大片免费视频hd| 午夜免费男女啪啪视频观看| 日本熟妇午夜| 成人三级黄色视频| 精品免费久久久久久久清纯| 欧美xxxx黑人xx丫x性爽| 伦精品一区二区三区| 日本免费一区二区三区高清不卡| 日日摸夜夜添夜夜爱| 免费人成在线观看视频色| 亚洲成人中文字幕在线播放| 国产视频首页在线观看| a级毛色黄片| 亚洲国产精品合色在线| 26uuu在线亚洲综合色| av卡一久久| 亚洲成人久久爱视频| 99精品在免费线老司机午夜| 又黄又爽又刺激的免费视频.| 国产成人精品久久久久久| 亚洲国产欧美人成| 国产视频内射| 国产精品国产高清国产av| 91aial.com中文字幕在线观看| 黄色视频,在线免费观看| 小说图片视频综合网站| 伦理电影大哥的女人| 女人十人毛片免费观看3o分钟| 精品久久久久久久久久免费视频| 精品不卡国产一区二区三区| 综合色av麻豆| 91久久精品国产一区二区成人| 精品99又大又爽又粗少妇毛片| 日本熟妇午夜| 国产私拍福利视频在线观看| 成人性生交大片免费视频hd| 好男人在线观看高清免费视频| 波多野结衣巨乳人妻| 国产成人freesex在线| 91午夜精品亚洲一区二区三区| 久久这里只有精品中国| 热99在线观看视频| 天美传媒精品一区二区| 干丝袜人妻中文字幕| 国产精品久久久久久久电影| 99热网站在线观看| 男女下面进入的视频免费午夜| 亚洲精品久久国产高清桃花| 国产成人精品一,二区 | 三级经典国产精品| 2021天堂中文幕一二区在线观| 国产精品,欧美在线| 久久午夜亚洲精品久久| 一级av片app| 久久精品91蜜桃| 在线a可以看的网站| 日韩中字成人| 亚洲av一区综合| 一级毛片久久久久久久久女| 99热这里只有是精品在线观看| 婷婷精品国产亚洲av| 国产一区二区在线观看日韩| 亚洲五月天丁香| a级毛片a级免费在线| 少妇人妻精品综合一区二区 | 一级毛片aaaaaa免费看小| av.在线天堂| 婷婷六月久久综合丁香| 精品一区二区三区视频在线| 成人美女网站在线观看视频| 九九久久精品国产亚洲av麻豆| 丰满乱子伦码专区| 成人鲁丝片一二三区免费| 国产熟女欧美一区二区| 国产高清激情床上av| а√天堂www在线а√下载| 久久精品夜色国产| 老女人水多毛片| 看免费成人av毛片| 波多野结衣高清作品| 天美传媒精品一区二区| 国产精品久久久久久久久免| 天堂网av新在线| 午夜a级毛片| 精品99又大又爽又粗少妇毛片| 久久久久国产网址| 长腿黑丝高跟| 最近视频中文字幕2019在线8| 老师上课跳d突然被开到最大视频| 老熟妇乱子伦视频在线观看| 99在线人妻在线中文字幕| 美女 人体艺术 gogo| 小蜜桃在线观看免费完整版高清| 成人一区二区视频在线观看| 男人舔女人下体高潮全视频| 91在线精品国自产拍蜜月| 日本撒尿小便嘘嘘汇集6| 亚洲国产日韩欧美精品在线观看| 久久婷婷人人爽人人干人人爱| 熟女人妻精品中文字幕| 成人综合一区亚洲| 亚洲国产精品sss在线观看| 亚洲精品粉嫩美女一区| 久久久精品大字幕| 伊人久久精品亚洲午夜| 国产午夜福利久久久久久| 少妇裸体淫交视频免费看高清| 晚上一个人看的免费电影| 国产成人a区在线观看| 又爽又黄a免费视频| 五月伊人婷婷丁香| 色视频www国产| 国产v大片淫在线免费观看| 国产三级中文精品| 一区二区三区四区激情视频 | 高清毛片免费观看视频网站| 亚洲国产日韩欧美精品在线观看| 亚洲自拍偷在线| 18禁在线无遮挡免费观看视频| 午夜福利视频1000在线观看| 亚洲在久久综合| 精品一区二区免费观看| 国产单亲对白刺激| 男女啪啪激烈高潮av片| 夜夜爽天天搞| 赤兔流量卡办理| 久久精品国产鲁丝片午夜精品| 色5月婷婷丁香| 国产精品一区二区三区四区免费观看| 成年女人看的毛片在线观看| 久久综合国产亚洲精品| 久久精品国产亚洲av天美| 午夜福利高清视频| 床上黄色一级片| 国产又黄又爽又无遮挡在线| eeuss影院久久| 国产v大片淫在线免费观看| 精品人妻偷拍中文字幕| 综合色av麻豆| 亚洲成a人片在线一区二区| 在线观看66精品国产| 国产三级中文精品| 男女视频在线观看网站免费| 亚洲欧美成人精品一区二区| 亚洲自偷自拍三级| 又爽又黄无遮挡网站| 边亲边吃奶的免费视频| 国产黄色视频一区二区在线观看 | 日韩视频在线欧美| 麻豆久久精品国产亚洲av| 国产国拍精品亚洲av在线观看| 国产美女午夜福利| 国产精品.久久久| 九色成人免费人妻av| 欧美另类亚洲清纯唯美| 成人毛片60女人毛片免费| or卡值多少钱| 联通29元200g的流量卡| 国产av在哪里看| 一个人免费在线观看电影| 亚洲天堂国产精品一区在线| 人人妻人人澡欧美一区二区| 欧美色视频一区免费| www.av在线官网国产| 欧美性感艳星| 亚洲av电影不卡..在线观看| 波多野结衣高清无吗| 国产精品人妻久久久影院| 亚洲精品乱码久久久v下载方式| 日本色播在线视频| 久久精品综合一区二区三区| 久久精品久久久久久久性| 亚洲av中文字字幕乱码综合| 国产午夜福利久久久久久| 国产精品爽爽va在线观看网站| 可以在线观看的亚洲视频| 国语自产精品视频在线第100页| 日韩成人av中文字幕在线观看| 看片在线看免费视频| 日本与韩国留学比较| 嫩草影院新地址| 亚洲人成网站在线播| 国产老妇女一区| 日韩精品有码人妻一区| 麻豆乱淫一区二区| 日韩,欧美,国产一区二区三区 | 能在线免费看毛片的网站| 中国美白少妇内射xxxbb| 少妇丰满av| 日本av手机在线免费观看| 亚洲国产日韩欧美精品在线观看|