Estimating the age of deciduous forests in northeast China with Enhanced Thematic Mapper Plus data acquired in different phenological seasons

Abstract This study investigated the ability of Landsat Enhanced Thematic Mapper Plus data acquired in leaf-on and leaf-off seasons to estimate stand age of Larix gmelinii and Betula platyphylla in northeast China. The relationships of six band reflectances, nine vegetation indices, and six texture measures with stand age were examined. Linear and multivariable regression models and multilayer perceptron neural network (MLP NN) were employed to estimate forest age based on these variables. The results indicate that reflectance in short-wave infrared bands and wetness are more significantly correlated with stand age in the leaf-on image, while reflectance in blue and green bands and greenness are more sensitive to stand age in leaf-off image. The MLP NN model can be effectively used to retrieve the stand age; the highest coefficient of determination and minimum root mean square error values of retrieved age are 0.47 and 21.3 years for Larix gmelinii, and 0.60 and 10.1 years for Betula platyphylla, respectively. The predicted age errors increased significantly when stand ages were > 100 and > 50 years for Larix gmelinii and Betula platyphylla, respectively. Remote sensing data acquired in the leaf-on season is more suitable for estimating forest age than that acquired in the leaf-off season over the study area.


Introduction
The carbon flux from the atmosphere to forests plays an important role in retarding the increase in atmospheric CO 2 concentration and climate change. 1 Forest age, which is related to time since disturbance or plantation, is a dominant determinant for the long-term trend of carbon exchange between forests and the atmosphere due to the age-related change in forest growth rate. 2 This phenomenon is caused by the changes in forest fundamental structure, which influences the ways of forest space occupied and carbon allocated. 3 Therefore, forest age is a crucial variable for quantifying carbon fluxes between the atmosphere and forest ecosystems. 4 It also plays an important role in estimating forest biomass with remote sensing data. 5 Forest age maps with high quality are urgently required for reliably quantifying regional and global forest carbon budget.
Forest age can be directly mapped according to the time of forest planted and disturbed using forest management inventory data and dendrochronology data. However, acquiring inventory data and tree core rings are often difficult, costly, and time-consuming, especially for natural forests. In addition, the data of forest age from inventories with broad age classes and low spatial resolution might lack the precision required for estimating carbon stocks and fluxes. For example, in many forest inventory datasets, forests are only grouped into different age ranges of 20 or 30 years and sometimes even just three succession stages (recent cut, secondary forest, primary forest) or five age classes (young, mid-aged, premature, mature, and postmature). 6 These broad age ranges and classes cannot satisfy the requirements for accurate calculation of forest carbon budget. Recently, remote sensing has been proved to be a useful tool for estimating forest stand age at reasonable cost, acceptable accuracy, and affordable efforts. The canopy structure, leaf area index, and background reflectance are main factors regulating the changes of remotely sensed reflectance with forest age. Other factors, such as sun and view angles and phenology, also have considerable effects on the reflectance and should be carefully addressed in the retrieval of forest stand age. 7 Forest age can be extracted using single or time-series remotely sensed data. Time-series remote sensing data are very effective for detecting forest disturbance history and forest age. However, the longest time-series of Landsat data can only be dated back to 1972. Forty-two years of remote sensing data are certainly not long enough for studying carbon budget and age-related growth decline of forest ecosystems. 8 In addition, the long time-series data with good quality is not available in many regions. Thus, we often should resort to the chronosequence method to retrieve forest age using a single remote sensing image. In this case, difficulties such as lacking suitable data and differences in illumination conditions and phenology could be avoided to some extent.
Remotely sensed spectral signals are not an explicit function of forest stand age. Rather, they are linked with structural factors related to age, such as biomass, leaf area index, density, basal area, and height. 9 Thus, it is possible for us to track forest succession stages and stand age according to changes in spectral signals. Several investigators have shown close relationships between forest succession stages and the spectral signals of Landsat Thematic Mapper (TM) data. 6,10,11 Jakubauskas and Price 9 successfully estimated the age (0 to 250 years) of lodgepole pine in the Yellowstone National Park using the data from visible, near-infrared (NIR), and middle-infrared bands of a TM image, with coefficient of determination (R 2 ) in the range from 0.62 to 0.90. Wulder et al. 6 used tasseled cap wetness indices calculated from a Landsat Enhanced Thematic Mapper Plus (ETMþ) imagery to estimate the stand age (4 to 20 years) of lodgepole pine in the Morice Forest District of British Columbia, Canada, with the root mean square error (RMSE) and R 2 of retrieved age 2.4 years and 0.68, respectively.
Texture measures are able to characterize vegetation structure and forest age. 11,12 Jakubauskas and Price 9 found that texture information from Landsat TM data was helpful to separate early forest succession sites in Yellowstone National Park. Lu 13 indicated that the combination of spectral and textural features of Landsat image improved the aboveground biomass estimation in Altamira, and texture was critical for mature forest biomass estimation. Texture measures were also proved useful for retrieving forest growing stock volume of oak in southern Liaoning, China, 14 and for delineating logged forests in Amazonia. 15 Most previous studies retrieved forest age using remote sensing data acquired in growing (leaf-on) season. Few studies tried to estimate the age of deciduous species using remote sensing data acquired in nongrowing (leaf-off) season. Forest basal area, which is related to the stand age, can be more effectively estimated using remote sensing images acquired in winter with snow ground cover than using those acquired in other seasons. 16 Snow produces a uniformly bright background that accentuates the signals of tree crowns and their shadows. Spatially integrated forest reflectance received by an optical satellite sensor during winter months is a complex mixture of sunlit/shaded tree canopies and forest floor, compounded by varying degrees of tree-totree self-shading and changes with the ratios of the four components. 16 In almost all wavelengths and all view angles, the reflected radiance of young deciduous stands is higher than the radiance of mature dense stands in winter because of more snow on the forest floor of young stands. 17 The differences in forest floor shadow or illumination are, in part, related to sun elevation angle, tree size, and stem density, and make it possible to estimate forest stand age using a single-date remote sensing image.
The relationships between forest age and spectral signals vary depending on the characteristics of the study area, forest types, and data acquisition time. No empirical models can be universally applicable for retrieving forest age from remote sensing data. China's forests account for a significant fraction of the global total and are playing an important role in global terrestrial carbon sink. 18 The carbon sink strength of forests in China significantly changes with forest age. 19 However, there is no national forest age map that can meet the need for reliably calculating the carbon budget of forests in China due to complex age structure of forests caused by intensive human and natural disturbances. Therefore, the development of algorithms to retrieve forest age based on remote sensing is of great scientific importance.
Studies have shown that reasonable accuracy of forest age retrieved from remote sensing data can be achieved using regression analysis, 6 traditional image classifications, 20 and machine learning algorithms. 21 New machine learning algorithms, such as artificial neural networks (ANNs), 22 support vector machines, 23 and random forest, 24 have emerged in the last years. These machine learning methods, without having to rely on statistical procedures or assumptions, have proven to be successful in land-cover classification, 25,26 but so far, few studies have compared their performance and usability in forest stand parameters retrieval. The multilayer perceptron neural network (MLP NN) is one of the most popular and successful ANNs architectures and has been successfully used for forest stand parameters estimation in recent years. 21,27,28 Thus, we chose MLP NN to test the ability of Landsat ETMþ images acquired in different seasons in indicating forest stand age.
In this study, we tried to explore the applicability of ETMþ remote sensing data acquired in leaf-on and leaf-off seasons in retrieving stand ages of two deciduous forest species (Larix gmelinii and Betula platyphylla), which are dominantly distributed within Daxinganling Area in northeast China. The main objectives of this study are (1) to test the ability of reflectance and calculated various spectral and texture indices to estimate forest age and (2) to assess the influence of acquisition time of remote sensing data on forest age estimation.

Study Area
The study area (123°19′ to 125°48′E and 52°09′ to 52°23′N) is a part of Daxinganling Forest Region, Heilongjiang Province, China (Fig. 1). Daxinganling region is well-known in China for extensively distributed forests. The elevation of study area ranges from 205 to 1399 m above the sea level. It has a cold temperate zone continental monsoon climate with an average annual rainfall of 350 to 500 mm and an annual mean temperature of −2.4°C. The winter here is harsh and long. Snow cover on the ground surface can last five months from November to March, with a depth of 30 to 50 cm (http://cdc.cma.gov.cn). The dominant tree species are Larix gmelinii and Betula platyphylla. The age classes (young, mid-aged, premature, mature, and postmature) are <40, 41 to 80, 81 to 100, 101 to 140, >141 years for Larix gmelinii and <30, 31 to 50, 51 to 60, 61 to 80, >80 years for Betula platyphylla, respectively. The understory plants are mainly Ericaceae and Cyperaceae. Larix gmelinii is normally in leaf at the beginning of May and defoliates at the end of September. Betula platyphylla begins to leaf out and defoliate in early May and mid-September, respectively. The growth of Larix gmelinii and Betula platyphylla normally ends during October and all leaves defoliate completely by the end of October. 29

Forest Resource Inventory Data
The forest resource inventory dataset collected in 1997 by local Forest Management Department was used in this study. This dataset includes detailed information of 29,722 subcompartment polygons, including stand age, percentage of dominant tree species, average diameter at breast height, average height, soil types, and elevation. They were recorded through field survey by local forest experts referring to disturbance history.
The forest resource inventory polygons satisfying the following conditions were used in this study: (1) forest stands are natural or conserving, (b) the ratios of dominant species (Larix gmelinii or Betula platyphylla) are >60%, (c) there is no evidence of recent damages caused by fire, insects, disease, or logging, (d) site class, a proxy for soil quality and forest productivity, is constrained to 3 and 4 to keep the polygons similar in environment conditions (site class in the range from 1 to 5 and mostly 3 and 4 in the study area), (e) the area of a polygon is >1 ha, and (f) the total area of polygons in a specific age is >100 ha. The ages of Larix gmelinii and Betula platyphylla range from 8 to170 years and 5 to 80 years ( Table 1). The ages of Larix gmelinii are almost normally distributed, peaking at ∼60 years. The Larix gmelinii forest ages were normally recorded at intervals of 5 years for planted forests and 10 years for natural forests and old planted forests by local forest experts. The numbers of polygons with the ages of 35, 45, 55, 65, 75, and 95 years are much smaller than those of polygons with ages of 30, 40, 50, 60, 70, 80, 90, and 100 years as natural forests are dominantly distributed in the study area (Fig. 2). The stands of Betula platyphylla are relatively young and the ages were normally recorded at intervals of five years. Above 90% of Betula platyphylla stands are at ages younger than 40 years (Fig. 2).  The images were geometrically corrected using a local topographic map as the reference. Topographic correction was conducted using a c-correction model. 30 The digital numbers were converted to surface reflectance in the following steps.

Image Preprocessing
First, Landsat ETMþ digital numbers (DN) were transformed into radiance using published postlaunch gains and offsets.
Second, the top-of-atmosphere (TOA) reflectance ρ TOA was calculated as 31 where L λ is the spectral radiance, d is the Earth-sun distance in astronomical unit, ESUN λ is the mean solar exoatmospheric irradiance, and θ s is the solar zenith angle in degrees. Third, the surface reflectance ρ s was retrieved as follows: 32 (3) where T g is the gaseous transmissivity due to the gases listed between parentheses, T RþA is the Rayleigh and aerosol transmissivity, ρ RþA is the Rayleigh and aerosols atmospheric intrinsic reflectance, and S RþA is the Rayleigh and aerosols spherical albedo. They were simulated using the 6 S radiative transfer code. 33 The two images (path/row 121/23 and 121/24) acquired on the same day were merged and then masked by the study area boundary. Snowpack with a depth of ∼13 cm existed in the image acquired on January 27, 2000. The influence of spatial change of snow property was not corrected here since it has minimal influences on spectral signals detected by the sensor compared to variations of shadow/illumination fractions and above-ground biomass. 16 Reflectance of bands 1 to 5 and band 7 and their different combinations were used here to estimate forest age, including brightness (B), greenness (G), and wetness (W), from the tasseled cap transform of ETMþ images, normalized difference vegetation index (NDVI), soil adjusted vegetation index (SAVI), normalized burn ratio (NBR), disturbance index (DI), integrated forest z-score index (IFZ), and moisture stress index (MSI) ( Table 2). These indices have been proved to be effective for mapping forest succession stage and forest age. [34][35][36] The reflectance of band 4 was used to calculate texture measures for estimating forest stand age because of its great variability with pixel brightness and high sensitivity to surface shadowing. 12  CON ¼ DIS ¼ HOM ¼ where pði; jÞ is the spatial co-occurrence matrix, i and j are the positional values within the matrix. The co-occurrence matrix was calculated in four directions (0, 45, 90, and 135 deg) for a 3 × 3 window with the distance 1 pixel and then summed to calculate each texture measure. 12,44 Means of reflectance, spectral indices, and texture measures for each polygon were calculated at two steps for analyzing their linkages with stand age. First, the spectral indices and texture measures were calculated for each pixel. Second, the means of reflectance, spectral indices, and texture measures of all pixels within a polygon were calculated.

Data Analysis and Model Construction
The relationships of reflectance, spectral indices, and texture measures with stand age were first analyzed by statistical approaches. The simple linear model (y ¼ a 0 þ a 1 x), logarithmic

Vegetation indices
Normalized burn ratio Normalized difference vegetation index Moisture stress index ETM5/ETM4 36 and 43 Note: ETM1, ETM2, ETM3, ETM4, ETM5, and ETM7 are the reflectance of bands 1, 2, 3, 4, 5, and 7, respectively. B, G, and W are the brightness, greenness, and wetness indices, respectively. B r , G r , and W r are the rescaled B, G, and W, respectively. B μ , G μ , and W μ are the means of B, G, and W, respectively. B σ , G σ , and W σ are the standard deviations of B, G, and W, respectively. b i and SD i are the mean and standard deviations of training pixels within an image for band i. b pi is the spectral value of pixel p in band i. NB is the number of spectral bands.
function (y ¼ b 0 lnðxÞ þ b 1 ), and quadratic model (y ¼ c 0 x 2 þ c 1 x þ c 2 ) were used to fit the relationships of stand age (x) with the variables calculated from remote sensing data (y). Then, the MLP NN model was employed to estimate forest stand age based on remote sensing variables. The usage of MLP NN includes two steps: training and prediction. During the first step, MLP NN utilizes a supervised error-learning technique called backpropagation for training the network. The learning process continues until the connection weights in the network have been adjusted so that the difference between output value and target value converges to an acceptable level. The second stage is to apply the trained network to new data to make a prediction. The polygon means of stand age and remote sensing variables were used as dependent and independent variables, respectively. To compare the performance of MLP NN, linear regression (LR) and stepwise multivariate regression (SMR) models were also used to estimate stand age. Over a large region, stand age is a continuous variable. However, in the forest resource inventory data used here, it is recorded as discretely incremental age classes, approximately at intervals of 5 or 10 years. Therefore, stand age retrieved from the forest inventory data are not really continuously distributed (Fig. 2). As an approximation, stand age was assumed as a continuous variable in the construction of models predicting ages here. Different combinations of reflectance, spectral indices, and texture measures were used as inputs into MLP NN, including (1) the reflectance of different ETMþ bands, (2) different spectral and structural vegetation indices, (3) different textural measures, and (4) all variables of reflectance, vegetation indices, and texture measures. A three-layer perceptron network was used here. The input was a vector of each combination of predictors and the output is the corresponding forest stand age. To evaluate more efficiently and accurately the performance of each model, we used a fivefold cross-validation technique, where the data were randomly partitioned into five subsamples: a single subsample was retained as the validation data for testing the model, and the remaining four subsamples were used as training data. The cross-validation process was then repeated five times, with each of the five subsamples used exactly once as the validation data. When testing was complete, we based our evaluation of each model on the predictions from all five testing samples.
In order to identify the importance of each variable for predicting forest stand age, sensitivity analysis (also known as independent variable importance analysis) was performed for all variables, which computes the importance and normalized importance value (NIV) of each predictor in determining the neural network based on the combined training and validation samples using the SPSS software. 45 The importance of an independent variable was calculated on the basis of how much model-predicted value changes if the value of this independent variable changes by AE10%, and the average of importance in the MLP NN model using the five subsamples was calculated. NIV was calculated as the ratio of the importance value of a specific independent variable to the largest importance value of all independent variables and expressed as percentage. The R 2 and RMSE of predicted forest stand age in comparison with those retrieved from the inventory data were computed to evaluate the overall accuracy of constructed models. The statistical analysis and model construction were performed with the SPSS statistics 17.0 46 and MATLAB®7.1 (The Math Works, Natick, Massachusetts) platforms.

Results
We tested the linear, quadratic, and logarithmic functions on forest stand age with each variable from leaf-on and leaf-off images, and the R 2 values and significance levels are shown in Table 3. The boxplots between forest age and remote sensing variables are shown in Figs. 3 to 5 (only the best and worst variables for each forest species and season are shown). For Larix gmelinii, the quadratic function is significantly superior to linear and logarithmic functions both in leaf-on and leaf-off images. For Betula platyphylla, the logarithmic function performed almost equally to quadratic function and better than the linear function in leaf-on and leaf-off images. Therefore, we selected quadratic function and logarithmic function for Larix gmelinii and Betula platyphylla, respectively.

Changes of Reflectance with Forest Stand Age
As for Larix gmelinii, in the image acquired on September 5, 1999, the reflectance of band 5 has the highest correlation with forest age (R 2 ¼ 0.138, p < 0.001), while the reflectance of band 2 is not significantly correlated with the forest age (p ¼ 0.402). The reflectance of bands 3, 5, and 7 shows similar trends changing with the age. It declines to minimums at stand age of ∼100 years and then increases. The reflectance of band 4 shows an opposite trend changing with stand age, peaking when stands are ∼34 years old. In the image acquired on January 27, 2000, the reflectance of all bands is significantly correlated with forest age, with R 2 in the range from 0.024 (band 7, p < 0.001) to 0.095 (band 1, p < 0.001). The spectral response of all bands to stand age is similar, with reflectance minimizing when stands were ∼90 to 100 years old. For Betula platyphylla, in the leaf-on image acquired on September 5, 1999, the reflectance of all bands decreases quickly with the age when forest is young and then changes marginally. The reflectance of bands 1, 2, 3, and 4 has stronger correlations with the forest age than that of bands 5 and 7 in the leaf-off image acquired on January 27, 2000. The highest R 2 values are 0.579 (band 5, p < 0.001) and 0.447 (band1, p < 0.001) for the leaf-on and leaf-off images, respectively.

Relationships of Vegetation Indices with Forest Stand Age
All spectral and structural indices are significantly correlated with the stand age except MSI in leaf-off image for Larix gmelinii (Table 3). In the leaf-on image, IFZ has the highest R 2 value (0.155) correlated with the stand age of Larix gmelinii, followed by W and B. Indices W, G, NBR, NDVI, and SAVI first increase and then decrease with the stand age of Larix gmelinii. The transition occurs at stand ages of ∼98 to 103 years. B, DI, IFZ, and MSI show opposite trends changing with stand age, decreasing first and then increasing. All indices with the exception of G, NDVI, and DI are more significantly correlated with stand age in the leaf-on image than in the leaf-off image. IFZ has the highest R 2 (0.145) correlated with stand age, followed by G, B, and NDVI in the leaf-off image. G and NDVI increase with stand age until stand age approaches 98 to 112 years and then decrease. In contrast, B, W, NBR, SAVI, DI, and IFZ decrease first and then increase with stand age.  In the leaf-on image, the R 2 values of W, NBR, NDVI, DI, and MSI against the stand age of Betula platyphylla are >0.3 (p < 0.001). G, W, NBR, NDVI, and SAVI show increasing trends with the stand age, while B, DI, IFZ, and MSI decrease with the stand age. The changes of W, NBR, SAVI, and MSI with the stand age of Betula platyphylla in the leaf-on image are opposite to those in the leaf-off image. In both leaf-on and leaf-off images, the correlations of most vegetation indices with the stand age of Larix gmelinii switch positively or negatively when the stand ages are in the range from 74 to 158 years. When stand ages are smaller than 80 years, the pattern might be different for each single vegetation index, but it is similar for Larix gmelinii and Betula platyphylla.    Figure 5 shows the changes of mean texture measures with the stand ages of Larix gmelinii and Betula platyphylla in two seasons. All textural measures are significantly correlated with the stand age except HOM in leaf-off image for Betula platyphylla (Table 3). In the leaf-on image, the changes of texture measures with forest stand age show less degree of parabola and logarithm than those of reflectance and vegetation indexes for Larix gmelinii and Betula platyphylla, respectively. The R 2 values of texture measures against the stand age are in the range from 0.046 (CON, p < 0.001) to 0.071 (ENT, p < 0.001) for Larix gmelinii and from 0.001 (HOM, p ¼ 0.066) to 0.230 (MN, p < 0.001) for Betula platyphylla in leaf-off image. As for Larix gmelinii, HOM slightly increases with the stand age until it approaches 100 years and then decreases with further increase of stand age. In contrast, other texture measures   show opposite trends changing with the stand age, decreasing first and then increasing. The transition occurs at stand ages ranging from 98 to 158 years. Texture measures exhibit similar changes with the stand ages of Betula platyphylla and Larix gmelinii before the stand ages approach the transition points. VAR, ENT, CON, and DIS perform better in indicating the stand age in the leaf-off season than in the leaf-on season, especially for Betula platyphylla.

Estimation of Forest Stand Age
The MLP NN method was used to develop models predicting forest age with different remotely sensed variables as inputs. Table 4 shows the statistics for the validation of age predicted using the MLP NN model with different multivariables as inputs. In the leaf-on image, the R 2 and RMSE values of predicted age are in the range from 0.13 to 0. 47    Larix gmelinii and from 0.18 to 0.60 and 14.5 to 10.1 years for Betula platyphylla, respectively. In the leaf-off image, the R 2 and RMSE values of predicted age are in the range from 0.11 to 0.26 and 27.5 to 25.1 years for Larix gmelinii and from 0.41 to 0.50 and 12.3 to 11.3 years for Betula platyphylla, respectively. Vegetation indices and reflectance perform better than texture measures in predicting forest stand age. The combination of all variables together can improve the prediction of forest stand age to some extent. The data acquired in the leaf-on season perform better in predicting stand age than the data acquired in leaf-off season. Figure 6 shows the comparison of stand age retrieved from the inventory data and predicted using the MLP NN model with all variables of reflectance, vegetation indices, and texture measures as inputs. The predicted age agrees well with the inventory age when the ages of Larix gmelinii and Betula platyphylla are in the ranges from 40 to 100 years and from 20 to 50 years, respectively. However, the model tends to overestimate the age of young stands and underestimate the stand age of old stands. The transition occurs when the ages are ∼70 years for Larix gmelinii and 30 years for Betula platyphylla. The RMSE of predicted stand age ranges from 12 to 62 years and 13 to 72 years for Larix gmelinii and from 7 to 38 years and 8 to 38 years for Betula platyphylla in the leaf-on and leaf-off images, respectively. As for Larix gmelinii, the RMSE of the predicted age is ∼20 years when the inventory age is <100 years and then increases dramatically. The RMSE of the predicted age of Betula platyphylla is ∼10 years when the inventory age is <50 years and then increases obviously.
Based on the sensitivity analysis for the multivariate MLP NN model, the most important predictors (NIV > 0.5) for the stand age of Larix gmelinii are ETM2, NBR, and ETM7 in the leaf-on image and IFZ, ETM1, and SAVI in the leaf-off image (Fig. 7). For Betula platyphylla, there are more variables that play relatively important roles in predicting the age (NIV > 0.5), including W, ETM5, ETM4, and NBR in the leaf-on image and G, ETM1, ETM7, ETM5, NBR, VAR, and DI in the leaf-off image.
Meanwhile, LR and SMR models were also used to estimate stand age of Larix gmelinii and Betula platyphylla for the five subsamples. The logarithm transformation was used to convert the forest age of Betula platyphylla to linear trend with predictors. For Larix gmelinii, the R 2 values of LR and SMR are 0.13 and 0.25 for leaf-on image, and 0.14 and 0.16 for leaf-off image, respectively. For Betula platyphylla, comparing with the LR model, the SMR model marginally improves the prediction accuracy of the stand age with R 2 increased from 0.52 to 0.54 in the leaf-on image and from 0.41 to 0.42 in the leaf-off image, respectively ( Table 4). The MLP NN model with all variables as input outperforms the LR and SMR models in predicting forest stand age.

Correlations of Reflectance, Vegetation Indices, and Texture Measures with Forest Stand Age
Understanding reflectance trajectories is important for detecting forest changes in satellite remote sensing applications and quantifying carbon sequestration by forest ecosystems. The nonlinear trajectories found in earlier studies in other geographic locations were most in leaf-on images. 27 This study revealed nonlinear relationships of remote sensing variables with forest stand ages of Larix gmelinii and Betula platyphylla in both leaf-on and leaf-off images, indicating the applicability of leaf-off images in estimating stand age. Forest spectral signatures were affected not only by the biochemical properties of leaves, but also by the spatial distribution and sizes of trees, which can be well sensed in the leaf-off images.
In this study, it was found out that reflectance of bands 1, 2, 3, and 4 shows much closer correlations with the stand age of Larix gmelinii and Betula platyphylla than that of bands 5 and 7 in the leaf-off image. In contrast, the reflectance of band 5 has much stronger linkages with forest age than the other bands in the leaf-on image. The correlations of reflectance in visible and NIR bands with forest age change with canopy density. 9,34 Shadow probably is an important factor altering the response of the reflectance in all bands to the changes in forest age. However, in the leaf-on image, the complexity of canopy and understory might impose noise on the relationship between forest age and the spectral signals. In the leaf-off season, the background is less complex and the shadow plays a dominant role in affecting remotely sensed signals under the condition of very low sun elevation, resulting in high sensitivity of the reflectance of visible and NIR bands to forest age. In the leaf-off image used here, the reflectance of visible and NIR bands is higher than that in the leaf-on image due to the high reflectance of snow on the floor and decreases with forest age until stand age approached >90 years. The reflectance of short-wave infrared bands 5 and 7 decreases with canopy and soil water content and canopy closure, and shows very low dynamic ranges under the low-illumination conditions with snow-covered ground in winter. Therefore, the reflectance of these two bands is more significantly correlated with forest age in the leaf-on image than in the leaf-off image.
In the leaf-on season, among three tasseled cap (TC) transform indices (B, G, W), W, which distinguishes among different classes of closed-canopy forests 10 and displays high correlation with stand attributes 47 and forest stand age, 11 outperforms B and G in depicting the stand age trajectories of Larix gmelinii and Betula platyphylla. Among three biophysical properties vegetation indices (NBR, NDVI, and SAVI), SAVI and NDVI, which are sensitive to living chlorophyll and moisture content of leaves and soils, show more sensitivity to stand age variations of Larix gmelinii and Betula platyphylla, respectively. Among three structure indices (DI, IFZ, and MSI), IFZ and DI, which were found very effective for mapping forest damage, show excellent ability to trace the forest stand age dynamics for Larix gmelinii and Betula platyphylla, respectively. The good predictors of stand age in the leaf-off image differ from those in the leaf-on image. The trajectories of W, NBR, SAVI, and MSI changing with stand age are opposite in leafon and leaf-off images for both Larix gmelinii and Betula platyphylla. This might be due to the sensitivity of the reflectance in near-infrared band (ETM4) to the snow background fraction in winter and insensitivity to forest floor conditions in the leaf-on season. 48 The reflectance, vegetation indices, and texture measures variables show better correlations with the age of Betula platyphylla than with that of Larix gmelinii. This mainly results from more complex tree structure of Larix gmelinii and its more interception of snow in comparison with Betula platyphylla in the leaf-off season. The observable fraction of sunlit snow dominantly determines the forest reflectance, especially for the visible spectrum region. 7 Texture measures indicate the heterogeneity in the tonal values of pixels within a defined area of an image. 12 For Larix gmelinii, the six texture measures in the leaf-on image show less obvious parabola trends changing with the stand age in comparison with reflectance and vegetation indices. The changes of texture measures with the age are similar for Larix gmelinii and Betula platyphylla in the leaf-on image. It should be noted that texture measures have complex features as they strongly depend on the window size as well as spatial and spectral resolutions. The appropriate window size for higher-resolution data might be more explanatory in regards to forest structure and may be more important than actual texture measure itself. 11 Here we used 3 × 3 window size in texture measures calculation for consistency in different forest species and images; however, more investigations should be conducted for investigating texture measures used to estimate forest stand age.
In the growing season, remotely sensed reflectance changes with phenology, resulting in corresponding changes in calculated vegetation indices and texture measures. The leaf-on image used in this study was acquired on September 5, 1999, just about two weeks prior to the initiation of defoliation. The color of Larix gmelinii and Betula platyphylla leaves in the remote sensing data acquisition time might have slightly changed. Therefore, it should be kept in mind that the relationships of reflectance, vegetation index, and texture measures with forest age identified here may not be directly applicable for remote sensing data acquired earlier in growing season.

Uncertainties in Estimated Forest Age
There are uncertainties in estimated forest age, which are possibly caused by the following reasons. The spectral signatures mixture significantly impacts the forest age prediction, especially for Larix gmelinii. Borders, roads, and dead trees within the forest stands may cause bias in the estimated age. The exclusion of stand edge pixels and resegmentation forest stand polygons can significantly improve the forest age estimation. 6,49 The low spatial resolution limits the accuracy improvement by combining texture and spectral information. 50 The wide age range and few polygons for young and old forest stands induce uncertainties in the MLP NN model. Of course, the coarse resolution and uncertainties of ages (at intervals of 5 and 10 years) in the inventory data might cause additional disagreement between predicted and observed ages. The changes in the accuracy of predicted forest age with the forest age were also reported for rubber forest in southern China. 49 Previous studies for Canada's boreal forests reported similar systematic biases of estimated age, i.e., an overestimation for young stands and an underestimation for old stands.
Young stands with open canopy and old stands with natural gaps are more easily impacted by the background (understory, leaf litter, grass, lichen, moss, rock, soil, snow, or their mixtures) due to the abundant understory and relative sparse canopy. The background reflectance inevitably impacts the relationships between forest stand attributes and satellite-measured reflectance, particularly in the low canopy cover of the tree layer. The difference in reflectance due to the understory vegetation can reach up to AE18% of the total reflectance in the red and up to 10% in the NIR regions of the spectra. 51 Recently, researchers have developed algorithms for extracting seasonal and spatial background reflectance values using the multiangle imaging spectroradiometer data and MODIS bidirectional reflectance distribution function model parameters product. 52,53 This would be helpful to explore the effect of background reflectance on these obtained relationships and should be further considered in forest stand attributes (e.g., forest age) retrieval.
According to our study, the leaf-off image covered with snow is not as good as the leaf-on image for predicting forest stand age. This might be caused by noise mentioned above, single winter remote sensing image used, and without consideration of spatial heterogeneity of snow cover. A previous study declared that the combination of two or more Landsat images acquired at different dates in winter might improve the estimation of forest basal area compared with the common approaches using either a single-date image or multitemporal summer images. 16 Also time series of winter satellite images combined with summer images would be more effective in estimating the forest succession stage after disturbance. 48 The present study shows the distinct trajectories of reflectance, vegetation indexes, and texture measures changing with the forest age, and the potential of mapping forest stand age over a large area with a snow-covered Landsat image. Remote sensing images with snow covered in the winter provide a unique view of the forest conditions. The application of multitemporal winter imagery and snow indices to explore the potential relationships between deciduous forest age and satellite data will be our next step.
The best variables for predicting forest stand age vary for different phenology and forest species. The important predictors in the MLP NN model are not identical to those identified in traditional statistical methods. The MLP NN provides dynamic NIVs as further data are fed to it. The better predictors identified in the univariate models and multivariate stepwise regression models do not certainly have higher NIVs. This is mainly owing to the ability of the MLP NN model to take into consideration the outliers and nonlinear interactions among variables and to reveal previously unrecognized and/or weak relationships between given input variables and an outcome. 45 Therefore, the MLP NN often includes variables that may not be significantly correlated with the predictand, as evidenced by the fact that the reflectance of ETM2 and ETM4 play important roles in the MLP NN model for estimating the ages for Larix gmelinii and Betula platyphylla in the leaf-on image, respectively, although the reflectance of these two bands is not so significantly correlated with the age. Validations showed that with reflectance, vegetation indices, and texture measures all together as inputs, the MLP NN model was able to improve the accuracy of forest age estimation for Larix gmelinii in both leaf-on and leaf-off images, and for Betula platyphylla in the leaf-off image (Table 4). This is largely due to different ability of most variables to indicate the stand age and different degree of interactions among them in leaf-on and leaf-off images ( Table 4, Figs. 3 to 5). For Larix gmelinii with a large age range and complicated canopy structure, reflectance, vegetation indices, and texture measures should be simultaneously applied for estimating forest age in both leafon and leaf-off images. As for Betula platyphylla with a small age range and less complicated canopy structure, forest age can be well estimated only with reflectance or vegetation indices in the leaf-on image.

Conclusions
In this study, we used Landsat ETMþ imageries acquired in different seasons to estimate the stand age of two typical forest species (coniferous Larix gmelinii and broad-leaved Betula platyphylla) in northeastern China. The ability of reflectance, vegetation indices, and texture measures to estimate forest age and the effects of remote sensing acquisition time on the retrieval of forest age were examined. Following conclusions can be drawn from this study: 1. There are nonlinear trajectories of remote sensing variables changing with the stand ages of Larix gmelinii and Betula platyphylla in both leaf-on and leaf-off images. These trajectories could be fitted by parabola and logarithm curves, respectively. 2. The correlations of remotely sensed variables with forest age depend on the data acquisition season. Reflectance in short-wave infrared bands 5 and 7, and wetness (W) have higher sensitivity to the variation of stand age in the leaf-on image, while the reflectance of visible and NIR bands 1, 2, 3, and 4, TC transform indices (B, G, and W), and NDVI are more significantly correlated with the forest age in the leaf-off image. Generally, texture measures are more sensitive to variation of stand age in the leaf-off image than in the leaf-on image, but do not show extraordinary sensitivity to the variation of forest stand age comparing with reflectance and vegetation indices. 3. The MLP NN model driven by remote sensing data is able to estimate forest age of young and middle-aged stands in the study area and outperform the LR and SMR models. As for old stands, the age predicted by the MLP NN model contains large uncertainties. The data acquired in the leaf-on season are more suitable for retrieving forest age than that acquired in the leaf-off season. Spectral properties overpass texture measures in predicting forest age in the study area.