Pines of Central Spain using QuickBird-2 Imagery and Classification and Regression Tree Analysis (CART)

Forest structural parameters such as quadratic mean diameter, basal area, and number of trees per unit area are important for the assessment of wood volume and biomass and represent key forest inventory attributes. Forest inventory information is required to support sustainable management, carbon accounting, and policy development activities. Digital image processing of remotely sensed imagery is increasingly utilized to assist traditional, more manual, methods in the estimation of forest structural attributes over extensive areas, also enabling evaluation of change over time. Empirical attribute estimation with remotely sensed data is frequently employed, yet with known limitations, especially over complex environments such as Mediterranean forests. In this study, the capacity of high spatial resolution (HSR) imagery and related techniques to model structural parameters at the stand level (n = 490) in Mediterranean pines in Central Spain is tested using data from the commercial satellite QuickBird-2. Spectral and spatial information derived from multispectral and panchromatic imagery (2.4 m and 0.68 m sided pixels, respectively) served to model structural parameters. Classification and Regression Tree Analysis (CART) was selected for the modeling of attributes. Accurate models were produced of quadratic mean diameter (QMD) (R 2 = 0.8; RMSE = 0.13 m) with an average error of 17% while basal area (BA) models produced an average error of 22% (RMSE = 5.79 m 2 /ha).


Introduction
Sustainable management of Mediterranean pine forests requires detailed and up-to-date information regarding structural parameters [1].Wood volume and biomass content in forest stands, calculated with structural indicators such as mean height and quadratic mean diameter, are basic data for administration of resources.Moreover, increasingly important and emerging environmental concerns related to habitat protection, carbon accounting, and biodiversity, make reliable knowledge of forest resources a requirement for national and international reporting [2].
In Spain, as in many other countries, accurate information of structural parameters is usually obtained via direct measurements by crews on the ground of systematically sampled field inventories, based upon a network of plots located on a regular grid [3] that is also subject to prior stratification.Field surveys are often costly and typically not spatially exhaustive.Field surveys are also often collected over a given re-measurement period, which can preclude adequate updating of information for periodic reports, and are of questionable validity over dynamic or non-merchantable forests.Despite these concerns, ground based inventories provide reliable and detailed information for development of models such as yield tables per species and given location.It is the difficulties in portraying these plot based measures spatially that for many applications limit the utility of this information to address more broad forest monitoring and reporting objectives [4], especially in heterogeneous forests.
Satellite imagery has been shown to support forest inventories of extensive areas by providing timely observation, increasing the accuracy of area estimates, producing wall-to-wall thematic maps, and providing inventory estimates with acceptable bias and precision [5].The spatially detailed information provided by high spatial resolution (HSR) imagery makes it an appropriate data source to aid in accurate estimation of structural parameters, and following suitable methods facilitates the characterization of subtle changes in forest structure through time [6].
The goal of this research is to explore the potential of HSR imagery to characterize forest structure in Mediterranean pines in the Central Range of Spain.Motivated by this purpose we examine the capacity of QuickBird-2 imagery to model the quadratic mean diameter, basal area, and number of trees per unit area at the stand level (as direct estimators of volume and biomass).Our specific objectives are:  To model the relation between structural parameters (quadratic mean diameter, basal area, and number of stems per hectare) measured via field sampling and a set of spectral and spatial variables derived from HSR multispectral and panchromatic imagery. To test and verify the ability of Classification and Regression Trees (CART) as the statistical technique for modeling structural parameters. To identify the image derived variables with the greatest informative capacity in the modeling of structural parameters, assessing in particular the inclusion of image textural metrics in the models.

Background
Space-borne optical remote sensing is a reliable source of information for assessment of forest characteristics over wide areas [7].The synoptic view and the regular acquisition cycle of image data, combined with the burgeoning selection of techniques available for attribute estimation, make remotely sensed data an appropriate and valuable source of data for assessment of forest condition and detection of change-offering information to augment costly and time consuming field campaigns for inventory update and re-measurement [8].

High Spatial Resolution (HSR) Imagery
Spatial resolution is an important consideration when using remote sensing for forest characterization [9].Currently the spatial resolution of systems frequently used for vegetation characterization range from coarse (e.g., 1 km of the Advanced Very High Resolution Radiometer) to very high (e.g., 0.4 m of the GeoEye-1 sensor).The adequacy of remotely sensed data for a specific purpose (e.g., attribute level: tree, stand, landscape, region) is conditioned by its spatial resolution, which is also inversely related to the extent covered by the image [10], also known as the image footprint.
Medium spatial resolution data with pixels sized 10-100 m (e.g., Landsat Thematic Mapper (30 m), ASTER (15 m)) are appropriate for characterization of forest condition [11] and monitoring of conditions and change at the forest stand level [12].Certainly a key to the applications and monitoring success of Landsat is the ability to capture conditions and dynamics that relate human interaction with terrestrial ecosystems.However, more detailed spatial data available since the launch of various commercial satellites (e.g., IKONOS in 1999, Orbview-3 in 2003) provide the opportunity for more precise depiction of forest parameters and are poised to reduce estimation errors of forest attributes to an acceptable level for operational applications [13].HSR imagery facilitates, for instance, the detection of individual tree characteristics [14], providing improved estimates of forest structural attributes [7].Panchromatic imagery, with fine spatial resolution (< 1 m) is particularly well suited for analysis of spatial relations through image texture measures [15,16].Texture measures enable the combination of spatial detail of panchromatic imagery with unique spectral information conferred by multispectral imagery serving to leverage complementary information [17] that can be employed separately or with a pan-sharpening approach [18,19].Spectral measures may be understood to inform on vegetation status, type, and condition with textural measures informing on vegetation structure.
Still, the dearth of established methods for image processing and the complex interactions between sun-sensor-surface geometry and forest structural characteristics [20], particularly in complex topographies, persist in making the use of HSR data challenging [6].HSR imagery acquired using space-borne platforms allows for data collection over remote areas, with predictable georadiometic qualities, and information content analogous to mid-scale aerial photography-commonly used for forest inventory purposes.Lidar (Light detection and ranging) technology has a demonstrated capacity to characterize forest structure [21][22][23][24] albeit with high costs persisting to limit operational, wide-area applications [25].Although lidar, with a capacity to collect highly detailed information regarding forest attributes, shows promise as a means to collect plot-like data for training attribute estimation algorithms applied to HSR imagery.

HSR Related to Forest Structure
The research literature is replete with studies relating forest structural parameters estimated from HSR satellite data (Table 1).Frequent techniques to obtain information from HSR images include crown isolation [26,27], shadow analysis [18,28], texture analysis [13,29,30], and geostatistical approaches [31][32][33].The capacity to characterize forest structural attributes typically decreases as crown closure increases [6], with an asymptotic relationship predictably emerging for vertically distributed attributes of forest structure [34].

Status in the Use of Remote Sensing for Estimation of Forest Structure in Spain
The Spanish Plan Nacional de Teledetección (PNT) is committed to acquiring complete national coverages of HSR satellite imagery annually [41] and to make data available for research at no cost.The acquisition phase started in 2008 [42], capitalizing upon archival data to backdate the database to 2005 coverage.Initial coverage consist of SPOT5-HRG XS+P (2.5 m) data, with other sensors being considered for future acquisitions [43].Access to this data represents a unique opportunity to incorporate HSR into Spanish forest inventories as an operational and low cost data source to meet a range of information needs.The data is to be collected with a primary focus on land-use land-cover change assessment [42], but capacity to generate information for forest monitoring and reporting can also be generated.
Encouraged by a readily available source of data there has recently been an increased interest by the Spanish research community in relation to remote sensing technologies and the potential application to forest environments, in particular the characterization of forest structure.Vázquez de la Cueva [44] explored relationships between forest structural attributes at the plot level (e.g., height, basal area, and crown canopy closure) and spectral information derived from Landsat Enhanced Thematic Mapper Plus (ETM+; 30 m pixel size) imagery combined with topographic data.The study considered three types of forest in Central Spain and applied a multivariate canonical ordination method.The author found a strong influence of vegetation type on the results, with a low percentage of variance explained precluding development of robust empirical models.Pascual et al. [45] used lidar data and a two stage object based methodology to characterize the structure of Pinus sylvestris L. stands in forests of Central Spain.Five structure types were defined based on height and density parameters.The median and standard deviation of height were found to be the most valuable for definition of structure types, with the approach developed being proposed for operational application suitable for inclusion in forest inventory procedures in support of forest management plans.Merino de Miguel et al. [33] investigated the strength of relations between dasometric parameters and textural variables in Pinus pinaster Ait.stands in Central Spain.The authors used geostatistical tools (i.e., variograms), calculated with orthophotography and IKONOS-2 imagery with original and degraded spatial resolutions.The authors found the strongest correlations when the variogram was calculated for spatial resolutions of 1 m and 2 m.As such, opportunities to further explore the capacity of HSR imagery to estimate a range of forest structural parameters remain.

Methods
Below, and in Figure 1, we summarize the approach implemented and the data utilized in this research.Forest structural attributes (QMD, BA, and N) are derived from data measured on the field through a process of geostatistical interpolation.Spectral and spatial variables from HSR imagery direct the delineation of stand-like areas for summarizing data.Statistical models linking forest parameters and imagery data are built with CART and validated with numerical and graphical tools.

Study Area and Field Data
The study focuses on pines in the Central Range of Spain (Figure 2), an area mainly dominated by P. sylvestris L., P. pinaster Ait., and P. nigra Arn.species.Two sites representing different forest conditions were chosen for availability of field data.Pinar de Valsaín (hereafter Valsaín) is a 7,627 ha forest of Pinus sylvestris L. on the North facing slopes of Sierra de Guadarrama (Segovia).It is a multifunctional forest (timber production, recreation, and protection) with an established management plan since 1889 that has evolved from a rigid to a more flexible scheme over the subsequent decades.
Management actions and recreational activities have had an impact on the forest structure [46].Valle de Iruelas (hereafter Iruelas) is a 5,483 ha forest of P. pinaster Ait., P. sylvestris L., and P. nigra Arn. in Sierra de Gredos (Ávila).It is also a multifunctional forest (wood, resin, and pasture production, recreation, and wildlife habitat).Although the first management plan was approved in 1886, historical circumstances prevented its implementation.The production of resin during the twentieth century favoured old growth development and a complex history of fires has also conditioned the forest structure.Systematic surveys based on ground sample plots are conducted periodically over the study sites measuring attributes including density, diameter at breast height (dbh), and height.For this study, data is from 2005 for Iruelas and 1999 for Valsaín, with the latter updated to 2004 conditions using a locally appropriate growth model following procedures recommended by the Spanish National Forest Inventory.The quadratic mean diameter (QMD) and basal area (BA) were calculated at each inventory plot (Equations (1-2)) where the total number of trees per unit area (N) was also available; expansion factors were used to scale values to a given area [47].BA and QMD are adequate attributes for volume modeling at the stand level.QMD was preferred over the arithmetic mean diameter as it has a stronger correlation to stand volume [48].Geostatistics provides a means for extrapolation of measured values to unmeasured points and areas, and facilitates the derivation of thematic layers for integration with other data [49].Kriging is a spatial interpolation method that yields the best possible estimation of the spatial variable of interest at every unmeasured point [50] and the error committed in the estimation is minimized and known at each point [51].In this study we mapped the forest variables of interest (QMD, BA, and N) measured in ground plots located over grids sided 150 m in Iruelas and 200 m in Valsaín into raster layers through a process of ordinary kriging.The relative standard error (i.e., the standard error of the kriged surface relative to the mean attribute value at the polygon level) was on average 15% for the QMD kriged layer and 25% for the BA and N layers, similar to the variability found for multiple plots found within the same polygon.More accurate averaging is facilitated, as sampling is complete and spatial correlation of plot values is accounted for.

HSR Imagery
QuickBird-2 is an Earth Observation satellite launched by Digital Globe in 2001, providing data in five spectral bands (Table 2).It has the capacity to be oriented and to capture images off nadir enabling a temporal revisit of 2-6 days depending on latitude [52].The pixel size of QuickBird-2 images is 2.4 m for the multispectral bands and 0.68 m for the panchromatic band (Table 2).
Two QuickBird-2 images, supplied in a georeferenced form by the data provider were used in this study, each covering one of the study sites (Figure 2, Table 2).Images were orthorectified with a Digital Elevation Model (DEM) derived from a contour vector map 1:10,000 (www.sitcyl.jcyl.es)and registered to aerial photography with 0.25 m pixels (www.sitcyl.jcyl.es).The multispectral and panchromatic bands were orthorectified separately with root mean square errors (RMSE) of 0.69-0.72 m (multispectral bands) and 0.66-0.81m (panchromatic band).Images were resampled with cubic convolution to 2.0 m (multispectral bands) and 0.6 m (panchromatic band) for alignment with the regionally appropriate coordinate grid (UTM 30N) and to facilitate integration with rasterized attributes.Atmospheric correction of the multispectral images was performed with the COST model [53] using water bodies as dark objects and the atmosphere-scattered path radiance L  p estimated with a relative spectral scattering DOS model ( −4 ) under very clear atmospheric conditions [54].

Image Segmentation
Image segmentation is the partitioning of images into uniform continuous spatial units [55].Through the application of automated algorithms the criteria for homogeneity can be defined by the user, based on parameters such as tone or spatial pattern.Image objects or segments composed of various pixels provide supplementary features for image analysis, not available in pixel based analysis, such as local statistical relations of digital numbers [37], shape, size or context.That is, once segments are produced, objects (i.e., trees or groups of trees) or spatially constrained summaries of the digital numbers within the segment may be used to provide representative segment-level information [39].In forest environments, the segments can often be considered as analogous to the manually delineated stands found in forest inventories [56].
Segmentation routines were applied to the QuickBird-2 images using Definiens Cognition Network Technology® [57,58].In the process of image segmentation the size of resulting objects is determined by the scale parameter and by the landscape characteristics; for instance a given scale value would produce larger objects in a homogeneous landscape and smaller objects in irregular areas.The scale parameter was 50 in Iruelas and 100 in Valsaín.Other settings guiding the segmentation routine include color-shape 0.8-0.2 and smoothness-compactness 0.5-0.5.The homogeneity criteria included the visible and NIR bands with similar weight, and an aspect layer derived from the DEM to incorporate topographic information as one of the possible structural driving factors [59] was weighted 0.1.

Image Texture Metrics
Image texture, defined by Haralick and Bryant [60] as "the pattern of spatial distributions of grey-tone", describes the relationship between elements of surface cover [61] and is one of the most valuable criteria in visual interpretation.The estimation of forest stand parameters is sometimes improved with a combination of spectral and spatial information [62] such as texture.Consequently a host of texture measures have been utilized to predict structural parameters in various environments [13,29,55,63,64] and has shown particular utility in complex structures such as tropical forests for above ground biomass estimation [17,65].We applied an approach for texture analysis based on measures derived from the Grey Level Coocurrence Matrix (GLCM) [66,67].The GLCM is a tabulation of how often different combinations of pixel grey levels occur in an image [68] at a specific distance and orientation (within a particular processing kernel, or analysis window).Texture analysis is a multiscale phenomenon [69] and choosing the right window size to capture meaningful local variance without generalizing unrelated features [13] is one of its key challenges [70].For selection of window sizes to calculate the GLCM texture measures we used the semivariogram approach [71,72].Semivariograms were calculated for image subsets over five experimental structural plots in Valsaín [73] and ten structurally different areas in Iruelas, identified with a combined approach based on inventory data and visual interpretation to cover all distinctive structural conditions.The range in the variogram indicates the distance beyond which pixel values are no longer correlated [71] and is an indication of the elements forming the texture present within the scene.The range is frequently associated with the most dominant elements in the scene, be it single tree crowns in open forests, or the canopy of groups of trees in close environments.Once the variograms were calculated, the range values were manually identified at the lag distance, where the variograms first flattened, corresponding with window sizes on the QuickBird-2 panchromatic band of 7 × 7, 9 × 9, and 13 × 13 pixels in Valsaín and 7 × 7, 13 × 13, and 23 × 23 pixels in Iruelas.We considered three GLCM texture variables, that is, Homogeneity, Contrast, and Entropy for each size of window (Small, Medium, and Large) (Table 3) based on their high values of correlation with structural parameters observed and pre-analysis investigations (results not shown).

Decision Tree
One option to identify relations between variables in multivariate data sets resulting from object analysis is the use of decision tree data analysis [37] also known as Classification and Regression Trees (CART).Regression trees identify relationships between a single continuous response (dependent variable) and multiple, continuous and/or discrete, explanatory (independent) variables, through a binary recursive partitioning process, where the data are split repeatedly into increasingly homogeneous groups (nodes), using combinations of variables (rules) that best distinguish the variation of the response variable.Tree models do not make assumptions regarding the distribution of the input data [74,75]; plus, they are able to capture non linear relationships between variables and are robust to errors in the input and results.Tree modeling is a nonparametric method which basic theory is reported in Breiman et al. [76].
CART approaches have frequently been used in the environmental remote sensing community for classification and mapping [77][78][79] for modeling [80][81][82] and for forest characterization [83].In the estimation of forest structural parameters with HSR satellite imagery, decision trees have been applied in diverse environments: Chubey et al. [37] used CART for analysis of percent species composition, crown closure, stand height, and age with IKONOS imagery based on analysis of objects in Alberta, Canada, obtaining the best estimations for species composition and crown closure.Goetz et al. [84] used IKONOS and shadow analysis to model and derive classified maps of canopy cover, with 97.3% overall accuracy, in Maryland, USA.Mora et al. [40] estimated mean height of forest stands in boreal coniferous forests in Yukon, Canada, obtaining a prediction accuracy of 53% and an RMSE of 2.84 m on stand height.All of the abovementioned approaches suggest local models for estimation of forest structural parameters as an alternative tool for alleviation of often costly and time consuming field inventories.

Applied Decision Tree
For development of decision tree models each segment was characterized with the mean and standard deviation of the reflectance and texture variables described above (Table 3), and the mean values of the kriged forest structural parameters (QMD, BA, N) and topographic orientation.These sets of data were input for the CART analysis in Matlab®.
Samples were randomly split into calibration (two thirds) and validation (one third) sets.The representativeness of the subsamples was tested with a Multi Response Permutation Procedure (MRPP) [85,86].This non-parametric method tests the hypothesis of no difference between two or more data sets for a range of parameters (i.e., the metrics used as inputs to the regression tree).To fit the model a cross validation process with ten iterations was performed; to avoid over-fitting we considered the establishment of a minimum number of cases in terminal nodes and pruning with the 1 SE rule [76].

Stand-Like Produced by Segmentation of the QuickBird-2 Imagery
Objects smaller than 0.5 ha produced in the process of segmentation were eliminated.Furthermore, screening outliers of reflectance and texture variables (i.e., segments which values were three or more standard deviations from the mean) enabled identification of objects that did not appear representative of known local forest conditions, typically corresponding with shepherding areas with buildings present in Valsaín and objects dominated by bare soil in Iruelas.Thirty nine such unusual objects were removed as outliers for subsequent analysis.Finally the number of objects preserved for modeling was 490, with an average area of 5.3 ha.Table 4 lists the statistical descriptors of the structural attributes (QMD, BA, N) and topographic parameter (aspect) at the stand-like level.Figure 3 illustrates the distribution of the structural parameters.Table 4. Statistical descriptors of structural (QMD, BA, N) and topographic (aspect) parameters of the stand-like objects obtained with the segmentation process and after removal of outliers.To fully capture the ecological meaning of the stand orientation and to avoid operational ambiguities we computed aspect values to be expressed as a non-polar complex number using the notation of Euler: Aspect = exp(−i × (θ − П/2)).

Regression Trees
Information regarding the calibration and validation subsamples is presented in Table 5.The MRPP test, performed including all stand level predictors, confirmed there were no significant differences between the calibration and validation datasets (p-value 0.77).Fitting all regression tree models was statistically significant (p-value < 0.001) and with high values of correlation (Table 6) between structural parameters and image predictors.To assess the performance of the models we applied them to the independent set of validation data, analyzing values of the Root Mean Square Error (RMSE) and correlation coefficient (R 2 ) (Table 6) and evaluating discrepancies between values measured on the field and values predicted by the regression tree models with the help of graphic tools (Figures 4 and 5).
Applied to the validation sample the models show varying strength of the relation between the structural parameters and the image variables used as predictors.The QMD model correlation value is the highest, followed by the BA model and with the N model ranking last (Table 6).The RMSE values, a means to measure the precision of the models, are moderate for QMD and BA, and relatively higher for N when a prediction of the exact number of trees is expected (Table 6).As practical decisions in forest management are often based on classes of attributes rather than exact values of structural parameters, we evaluated the performance of the CART model to classify values of N. The measured number of trees per unit area (N) was classified into density categories ranging from open (N < 150) to closed (N > 500) categories.The CART model classified 70% of the stand-like segments in the correct group, with all other segments classified in an adjacent class.The average relative error of the models was also evaluated as the percentage of RMSE respect to the average measured parameter (Table 6).Scatter plots in Figure 4 illustrate the relation of observed values of QMD (a), BA (b) and N (c) versus the corresponding estimated values of the validation subsample (n = 163).The QMD model performs with very good accuracy for the smaller diameters, with points close to the 1:1 line, and more randomly spread to both sides for larger diameters.The BA model depicts a similar but less accurate pattern, while the N model shows increasing disagreement of observed to modeled values at the more dense stands.Noteworthy is a tendency of underestimation for parameters at high values (QMD ≥ 1.2, BA ≥ 50, and N ≥ 600), likely as an expression of the well known saturation of optical sensors at increasingly high biophysical parameter values [34,87].This kind of error is important to note with reference to volume and biomass estimation, since larger trees contribute more to these estimates [88], but it is of minor importance in this particular area where few stands are over the thresholds mentioned above (Figure 3; Table 4).A closer look at the residuals confirms the relative precision of the QMD model (Figure 5(a)); an assessment of relative errors revealed that the relative error committed is below 20% in 76% of the validation sample (n = 123).A comparison of 5 cm diametric classes between the estimated and observed data indicated an agreement in 53% of the stand-like segments, with 19% falling in the adjacent class.Furthermore, the random distribution of residuals in the most frequent classes (0.30-0.40) leads to an almost complete compensation of the average error.This optimistic result should be carefully considered, as averaged values over areas of different sizes could lead to miscalculations.The residuals in the BA model look randomly distributed (Figure 5(b)), but there is a higher number of underestimates (57% of the validation sample) and in these cases the absolute value of residuals is higher.In the N model 55% of the validation segments are underestimated; a tendency to underestimate lower values and overestimate higher densities is observed.To reduce over-fitting and to make the models practical and operationally viable we established a minimum number of cases in terminal nodes (n = 80).Furthermore, examining the terminal nodes average values and the improvement of intra-group variance they represent from father nodes (i.e., decreased variance) appropriate pruning levels were determined.With these premises the number of terminal nodes obtained was between seven (for the QMD and BA models) and eight (for the N model) (Figure 6; Table 7).
The most relevant predicting variables determining decisions in the regression tree models are shown in Table 7. Noteworthy is the primacy of stdev B1 (standard deviation of blue reflectance) which enters all models in first place.All other reflectance bands (green, red and near-infrared) did also determine some branch rules (Figure 6).Among textural variables, contrast and entropy of various window sizes were the more relevant; homogeneity was not included in decision rules.A total of five or six variables were included in each of the models.

Discussion
Structural parameters such as quadratic mean diameter, basal area, and number of trees per unit area of Mediterranean pines in Central Spain have been modeled with regression trees and with HSR reflectance and texture metrics from QuickBird-2 imagery as model inputs.Results, although limited by uncertainties in the reference data and processing techniques, show reasonable accuracy (R 2 = 0.8) and precision (estimation relative error ~17%) for the QMD model and robust models (R > 0.7) for BA and N but with higher estimation relative error (22-31%).
Management plans were initiated in Spanish forests more than a hundred years ago [89].Albeit the early start, only 19% of the treed forest area in Spain is currently governed by a management plan under formal implementation [90].Often noted as a primary reason for this unfavorable proportion, is the high cost of field inventories, limiting surveys to forests with potential to produce economic revenue.However, with the increasing concern over environmental issues, current forest inventories are aimed at informing a variety of long-term objectives including biodiversity, carbon accounting, habitat protection and sustainable timber production [91].Remote sensing can contribute to the ability to produce timely, cost efficient inventory estimates via image segmentation for stand delineation [45] and statistical modeling for assessment of attributes with acceptable precision [5].HSR satellite sensors emerged a few years ago as promising data sources for forest inventory [6,92] providing consistent and frequent imagery.Our study demonstrates that in Mediterranean pines of Spain QuickBird-2 imagery and CART modeling would be useful and affordable for assisting in the assessment of forest areas with a variety of objectives (e.g., recreation, carbon storage), though caution is required to deal with inherent modeling uncertainties.Although remote sensing is not expected to replace completely field measurement any time in the near future [5] it would facilitate planning and management with realistic goals.
Among the strengths of HSR imagery is the high geometric fidelity [93] and the possibility of identification of individual elements such as trees or groups of trees.The unique capabilities of the QuickBird-2 instrument are exploited here by including texture metrics in the modeling, as image texture is influenced by biophysical parameters like crown diameter, distance between trees, tree positioning, LAI, and tree height.The historic limited use of texture parameters is often indicated as related to a paucity of appropriate software tools [94] and is being progressively overcome.Alternately, for monitoring programs with various dates of imagery and more than one scene, off-nadir view angles and differing solar and atmospheric conditions should be considered [20] as they may pose analysis difficulties.
Heterogeneous environments typically require a dense network of sample plots for an adequate assessment of varying conditions [95]; likewise, the capacity of a grid of inventory plots to capture the diversity of Mediterranean forests could be argued.With the complete coverage offered by remotely sensed data, selective sampling may become unnecessary, for instance if imputation techniques are applied.Furthermore, in applications where sampling is needed, segmentation of HSR images helps the design of sampling units by automatically and consistently defining homogeneous areas [96], otherwise delineated with human expert and costly effort.If adequately trained, segmentation algorithms have the ability to semi-automatically divide images into structurally homogeneous areas only requiring human revision [25], that can be used as strata to optimize the field sampling design [97] and also allowing the reduction of sample collection needs.
Tree models are easily interpreted and applied, with few statistical requirements imposed that make it an appropriate method of estimation in forest environments.Employing data from managed stands' field inventories in the support of modeling efforts has an intrinsic limitation related to the dearth of measurements of small trees; this circumstance is possibly related to a bias of the data considered as truth, and could partly excuse the underestimating trend of our models.All sources of uncertainty should be thoroughly considered for aiding the interpretation of modeling results.Our calibration dataset consisting of 327 stands is relatively large (66% of the sample) as the accuracy of decision tree models tends to increase with increasing calibration sample size [70].Mora et al. [40] in Yukon (Canada) demonstrated that a smaller calibration dataset (30% of the sample) could perform adequately if there were difficulties to obtain reference information, making this method an even more appealing tool for inventory.With a simple structure, that is, low number of rules and final nodes, CART constitute a practical and parsimonious tool for classification of stands for management or planning.The acquisition of periodic HSR coverage of the whole territory by the PNT poses an unprecedented opportunity to use remote sensing for assessment of the structure of Spanish forests that managers should strongly consider.

Conclusions
High spatial resolution (HSR) satellite imagery, such as QuickBird-2, has information content enabling the modeling of structural parameters for the pine forests of Central Spain.In this research the quadratic mean diameter (QMD), basal area (BA), and number of trees per hectare (N) of pines in the Central Range of Spain were modeled at the stand level with classification and regression trees (CART).Models were produced with average estimation errors suitable for planning purposes: predictions of QMD had an average error of 17% and BA an average error of 22%, while N was correctly classified in 70% of the cases.Although some refinement of the techniques applied here is possible to support operational activities, this study has demonstrated that following the selection of appropriate statistical tools combined with the periodic acquisition of HSR imagery by the Spanish Plan Nacional de Teledetección (PNT) could be of great value to the forest community as a low cost option to support planning activities.Additional stakeholders could also be accommodated and supplied with wide-area estimates of forest structural attributes following the methods suggested in this research.The capacity to revise the estimates with new plot data in subsequent years and to incorporate depletions using change detection procedures also points to additional utility and value that can be created from the national PNT image collections.

Figure 1 .
Figure 1.Schematic methodology followed in the study.

Figure 2 .
Figure 2. Location of the study sites.Insets show QMD values as kriged from inventory plots in the treed areas of Valle de Iruelas and Pinar de Valsaín.Subset areas covered by 834 plots in Valsaín and 661 plots in Iruelas were investigated in the study.

Figure 3 .
Figure 3. Distribution of the structural parameters (QMD, BA, N) in the stand-like polygons produced with the segmentation of the satellite images.Note that QMD graph bins are not all equal.

Figure 4 .
Figure 4. Plot of the observed structural parameters QMD (a), BA (b), and N (c), versus estimated values for the validation subsample (n = 163).

Figure 6 .
Figure 6.Example of a regression tree model of QMD.Hollow boxes represent branch rules; elements fulfilling the rule go to the left, the rest go to the right.Values of terminal nodes average QMD of elements in the group.

Table 2 .
Characteristics of the satellite imagery used in the study.

Table 3 .
Attributes used for modeling.The mean and standard deviation of each of these attributes was de facto used in the decision trees.

Table 5 .
Number of samples used for calibration and validation of the CART models.

Table 6 .
Fitting and performance results of the regression tree models for QMD, BA, and N.

Table 7 .
Relevant predictors in regression trees of QMD, BA and N and number of terminal nodes.