Classification of forest land attributes using multi-source remotely sensed data

International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

Contents lists available at ScienceDirect

International Journal of Applied Earth Observation and Geoinformation journal homepage: www.elsevier.com/locate/jag

Classification of forest land attributes using multi-source remotely sensed data Inka Pippuri a,∗ , Aki Suvanto b , Matti Maltamo a , Kari T. Korhonen c , Juho Pitkänen c , Petteri Packalen a a

School of Forest Sciences, University of Eastern Finland, P.O. Box 111, 80101 Joensuu, Finland Blom Kartta Oy, Kauppakatu 15, 80100 Joensuu, Finland c Natural Resources Institute Finland, P.O. Box 68, 80101 Joensuu, Finland b

a r t i c l e

i n f o

Article history: Received 26 March 2015 Received in revised form 2 July 2015 Accepted 3 July 2015 Keywords: Classification Forest land Landsat LiDAR Site type Surface model

a b s t r a c t The aim of the study was to (1) examine the classification of forest land using airborne laser scanning (ALS) data, satellite images and sample plots of the Finnish National Forest Inventory (NFI) as training data and to (2) identify best performing metrics for classifying forest land attributes. Six different schemes of forest land classification were studied: land use/land cover (LU/LC) classification using both national classes and FAO (Food and Agricultural Organization of the United Nations) classes, main type, site type, peat land type and drainage status. Special interest was to test different ALS-based surface metrics in classification of forest land attributes. Field data consisted of 828 NFI plots collected in 2008–2012 in southern Finland and remotely sensed data was from summer 2010. Multinomial logistic regression was used as the classification method. Classification of LU/LC classes were highly accurate (kappa-values 0.90 and 0.91) but also the classification of site type, peat land type and drainage status succeeded moderately well (kappa-values 0.51, 0.69 and 0.52). ALS-based surface metrics were found to be the most important predictor variables in classification of LU/LC class, main type and drainage status. In best classification models of forest site types both spectral metrics from satellite data and point cloud metrics from ALS were used. In turn, in the classification of peat land types ALS point cloud metrics played the most important role. Results indicated that the prediction of site type and forest land category could be incorporated into stand level forest management inventory system in Finland. © 2015 Elsevier B.V. All rights reserved.

1. Introduction Forest land and site type classifications are used at national and international levels to monitor the amount, quality and state of forests. In forestry, different classifications of forest land attributes are also important, as they form the basis for silvicultural operations. In Finland, forest land and site type information is collected in the National Forest Inventory (NFI) and in regional stand-wise forest inventories made for forest planning. The traditional stand-wise inventory, based mainly on field work, has been recently replaced with an inventory method using airborne laser scanning (ALS), aerial images and field measurements (Maltamo and Packalen, 2014). ALS-based inventory is accurate, less time-consuming and less expensive than the traditional stand-wise field methods (e.g. Næsset, 2004; Packalén and

∗ Corresponding author. E-mail address: inka.pippuri@uef.fi (I. Pippuri). http://dx.doi.org/10.1016/j.jag.2015.07.002 0303-2434/© 2015 Elsevier B.V. All rights reserved.

Maltamo, 2007). In the area-based ALS approach metrics calculated from height, density and intensity distribution of ALS points are used to predict stand attributes at plot or substand level (Holmgren, 2004; Næsset, 2002). Even though the stand attributes (volume, basal area, mean height and mean diameter of growing stock) can be predicted with high accuracy, classification of forest land and site types has only been tested in a few studies (Holopainen et al., 2010; Vehmas et al., 2011). ALS-based surface models have been found to be useful in many forestry applications. Despite this the potential of ALS-based surface metrics have not been extensively studied and e.g. in practical ALS-based forests inventories only ALS point cloud metrics are used. Surface metrics can be calculated, for example, from ALS-based terrain (DTM), surface (DSM) or canopy height (CHM) models. ALS-based DSM or CHM have been used in individual tree detection (Hyyppä and Inkinen, 1999), stand delineation (Koch et al., 2009) and evaluation of canopy cover and leaf area index (Korhonen et al., 2011), but also, e.g., in classification of forest types and estimation of forest attributes (Van Aardt et al., 2008,

12

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

change detection (Vastaranta et al., 2012) or prediction of spatial pattern of trees and need for forest management operations (Pippuri et al., 2012). ALS-based DTM metrics have only been utilized in a few forest studies. One example is Racine et al. (2014), who predicted stand age using ALS-based DTM metrics (e.g. elevation, slope, aspect, catchment area, solar radiation and wetness index) and point cloud metrics as structure metrics. Korpela et al. (2009) also used ALS-based surface metrics in the classification of mire site types. Metrics calculated from different surfaces may also describe forest land attributes better than straightforward ALS point cloud metrics. Most remote sensing-based forest classification studies have focused on the utilization of passive satellite images. There are also studies where forest land and different forest types have been discriminated using ALS data and some of these have utilized ALSbased surfaces, such as DTM, DSM, CHM and intensity rasters in classification (e.g. Antonarakis et al., 2008; Brennan and Webster, 2006; Charaniya et al., 2004). In some land use studies wetlands and swamp forest have also been mapped using remotely sensed data (e.g. Maxa and Bolstad, 2009; Sader et al., 1995; Townsend and Walsh, 2001). As far as we are aware, there are no ALS-based studies where mineral soils and peat lands have been separated from each other, but Dirksen (2013) discriminated swamp forests and non-paludified forests using ALS data. ALS data has only been used to classify forest site type (understory vegetation-based) in a couple of studies. Vehmas et al. (2011) classified forest site types in mineral soils using ALS point cloud metrics and Korpela et al. (2009) classified, e.g., mire site types, mire habitats and nutrient status based on point cloud metrics and ALS-based DTM metrics. Recently, the potential of ALS data to identify herb-rich forests has also been studied (Vehmas et al., 2009). Worldwide, a more common method for site type classification is site index, which is based on the dominant height of growing stock at a certain age. For instance, Gatziolis (2007) estimated dominant height and site index using an individual-tree-based ALS method and Packalén et al. (2011) using an area-based method. Forest land and site type information collected in NFI could be used as training data in practical ALS-based forest inventory. Earlier Hollaus et al. (2007), Maltamo et al. (2009) and Tuominen et al. (2014) successfully combined NFI plot data and ALS data for the estimation of growing stock attributes. Using comprehensive field data from NFI could increase the cost-effectiveness of ALS-based inventory, since separate field data would not need to be collected. In addition, the accuracy of the inventory could improve when forest land and site type classification is done by an NFI field expert. In this study the following classification schemes are studied: land use/land cover (LU/LC) classification using both national (Finnish) classes and United Nations Food and Agriculture Organization (FAO) classes, main type, site type, peat land type and drainage status. In the NFI of Finland both national and FAO LU/LC classification is used to assess the amount of forest area. National LU/LC classification is based on the purpose of land use (productive forest, poorly productive forest, open land, other forestry land, agricultural land, built-up land, etc.) and site productivity. Thus, the national LU/LC classification system can be regarded as a two-phase system: first classification by land use in forestry land, agricultural land, and built-up land classes, and then by site productivity where the forestry land is divided into productive forest, poorly productive forest, open land, and other forestry land sub-categories. The FAO LU/LC classification, in turn, is based on land use and crown cover and height of trees in their maturity stage. FAO LU/LC classification is used in the assessment of global forest recourses (FAO, 2006). In Finland forest land is classified into mineral soil or peat land (main site type) based on the composition of an organic layer of soil (peat or not) or proportion of peat land vegetation. Forest land

is also classified into uniform forest site types based on site fertility and wood-producing capacity. Site classification uses Cajander’s (1926) forest site type theory, which is based on understory vegetation. Mires constitute 28% of land cover in Finland (Korhonen et al., 2008) and are important habitats from the point of view of both forestry and biodiversity. They can be discriminated into dozens of peat land types (Laine and Vasander, 2008), but more general categorisation divides them into spruce- and pine-dominated and open peat lands. Peat lands are also often classified into corresponding site type classes with mineral soils. Forest land can also be classified based on drainage status, which describes the draining state of a forest stand. Draining status can be divided into undrained and drained (by ditches) mineral soils and peat lands. The main aim of this study was to (1) examine the classification of forest land attributes using ALS-based data, satellite images and NFI plots as training data and to (2) identify best performing metrics for classifying forest land attributes. It was of special interest to test different ALS-based surface metrics in the classification of forest land attributes. The predicting power of common ALS point cloud metrics and spectral metrics of satellite images were also evaluated.

2. Materials and methods 2.1. Study area and field data The study area (∼6500 km2 ) is located in Southern Finland (Fig. 1) and consists mainly forests, but also some lakes, fields and urban areas. Forests can be considered as typical managed boreal forests in Finland and are mainly dominated by coniferous species: Norway spruce (Picea abies (L.) Karst.) and Scots pine (Pinus sylvestris L.). Deciduous species are usually found as a minority in mixed forests. The study area includes both mineral soils and peat lands and forests with different site types and developing stages. Field data of 828 plots in 164 clusters (9–12 plots per cluster) were collected in the 10th and 11th cycles of the National Forest Inventory during 2008–2012. In NFI, a forest stand is defined as a homogeneous parcel regarding site and growing stock variables and it should usually be at least 0.25 ha in size. Land use and site variables of forests are always collected for the forest stand, or the parcel of land use class, in which the centre point of the plot is located. Additional forest stands are described if any trees in additional stands are selected to be measured. In Southern Finland, tree selection was done by a Bitterlich relascope with a basal area factor of 2 and a maximum radius of 12.52 m for the plot of trees. The proportion of the centre point’s stand within a circle with this maximum radius was recorded. In this study, we used plots that had only one forest stand or national LU/LC class within the 12.52 m radius. Plots in lakes/rivers and plots where the height of vegetation was clearly different compared to laser-based canopy height were excluded (caused, e.g., by clear cuts between the field and ALS data collections). Reason for this was that we wanted to focus on classification of forest land attributes (waters are very easy to detect using satellite data) and changes in forest between field inventory and remote sensing data acquisition affects remote sensing metrics in an apparent manner. Therefore, it does not make sense to include sample that are clearly wrong. Field data collected in NFI of Finland by experienced field workers can be expected to be very accurate evaluation of forest land attributes. The following classifications determined for each plot were used: national and FAO LU/LC class, main type, site type, peat land type and drainage status. Some of the classes only had a small number of observations, and therefore a few classes were merged. Despite this, the number of observations between different classes varied a lot. In the final classification schemes the national LU/LC classification included forestry land, agricultural

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

13

Fig. 1. Location of the study area.

land and built classes. Forestry land included productive forest land (98% of forestry land observations), poorly productive forest land, open land and other forestry land (e.g. forestry roads and timber storage areas). Class built included built-up land, roads

and power supply lines. The FAO LU/LC classes were merged into two classes, i.e. forest and non-forest, the latter including other wooded land and other land. The main type consisted of mineral soil and peat land classes. Forest site type classes were rich

14

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

Table 1 Classification schemes and included classes.

2.2. Reference plot types

Classification

Classes

National LU/LC FAO LU/LC Main type Site type Peat land type Drainage status

Forestry land, agricultural land, built Forest, non-forest Mineral soil, peat land Rich, medium, poor Spruce, pine, open peat land Drained, undrained

(as Oxalis–Maianthemum and Oxalis–Myrtillus types), medium (Myrtillus type) and poor (Vaccinum, Calluna and Cladina types and rock and sand sites) (Cajander, 1926). Each site type class indicates the corresponding fertility on mineral soils and peat lands. Peat land classes were spruce peat land, pine peat land and open peat land. Drainage status was divided into undrained and drained classes. The drained class included drained (ditched) mineral soils and recently drained, transforming and transformed peat lands. Final classification schemes are shown in Table 1 and the number of plots or segments in different classes in Fig. 2.

Two different reference plot types were used in modelling. First we used fixed 12.52 m-radius circular plots (the same as the maximum size for tree measurement in NFI). The centres of the plots were located in the field using GPS (global positioning system) positioning according the instructions of NFI of Finland. Another test was to use segments different in size. The segments were created by generating a 50 m buffer around the plot centre and delineating a homogeneous segment within the buffer zone (Fig. 3). Homogeneous segments were manually delineated/digitized around the plot based on canopy height model. The smaller circular plot with the 12.52 m radius was always included in the segment (because we rejected divided plots from the data), hence both plots and segments have the same stand information (categories) from the field assessment. We assumed that the larger segment represents better the forest stand or land use area described in the field measurements compared to a small circular plot. This is because forest land attributes are defined from the forest stand or parcel of land use class (not from the small plot), in which the centre point of the plot is located. Use of larger segments also reduces the random variation

Fig. 2. Number of plots or segments in different classes in classification of national and FAO LU/LC class, main type, peat land type, site type and drainage status.

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

15

calibration was carried out by taking a systematic sample of pixels from the overlapping area of the images and fitting robust regression models band by band (Tokola et al., 1999).

2.4. Predictor variables

Fig. 3. Reference plot types: circular plot with 12.52 m radius, segment and 50 m buffer, in which the segment was delineated based on canopy height model.

inside segments and might better detect e.g. sparse ditch networks compared to smaller circular plots. 2.3. Remotely sensed data Leaf-on ALS data were collected on 27 June to 4 July 2010 using an Optech ALTM laser scanning system at a flight altitude of 2000 m above ground level and with a maximum scan angle of 30◦ . This resulted in a swath width of 1050 m and a nominal sampling density of about 0.54 measurements per m2 . The Optech ALTM Gemini laser scanner captures one to four range measurements for each pulse, but the measurements were classified to represent first, last and intermediate echoes. The first echo data contained the echo categories ‘first of many’ and ‘only’, while the last echo data contained ‘last of many’ and ‘only’ echoes. A digital terrain model (DTM) was generated from the ALS data. First, the laser points were classified to ground points and other points (method explained by Axelsson (2000) by National Land Survey (NLS) of Finland and then a DTM raster was created from the ground points by taking the mean height of the points within each raster cell. The orthometric heights of laser hits (z value) were converted to above-ground heights by subtracting the DTM at the corresponding location. A canopy height model (CHM) was generated from the first echo data by taking the maximum height at the above-ground scale within a 1.6 m radius from the centre of a pixel. A digital surface model (DSM) was generated by summing up the pixel values of the DTM and CHM rasters. The pixel size of all ALS surface models was 1 m. Satellite data consisted of two Landsat 5 TM images acquired on 28 June and 7 July 2010. The spatial resolution of the images was 30 m. Landsat 5 TM imagery has seven bands: band 1 (blue) 0.45–0.52, band 2 (green) 0.52–0.60, band 3 (red) 0.60–0.69, band 4 (near infrared) 0.77–0.90, band 5 (shortwave infrared) 1.55–1.75, band 6 (thermal infrared) 10.40–12.50, and band 7 (reflective infrared) 2.08–2.35. There were some clouds in another image which were masked out with the Fmask software (Zhu and Woodcock, 2012). Images were then calibrated to another image to make the images radiometrically more similar. The radiometric

The metrics used as predictor variables were calculated for the area of reference plots and segments from ALS point cloud data, ALS surface models and Landsat images (Table 2). The principles of the area-based method (Næsset, 2002) were used in the calculation of ALS point cloud metrics. First, the height distributions were calculated for each plot and segment using the above-ground heights of the echoes. Distributions were calculated separately for first (f) and last (l) echoes and all laser hits were included as well as ground hits. Weighted height percentiles of 1, 5, 10, 20, . . ., 90, 95, 99 (h1, . . ., h99) were computed by summing the heights at above-ground level. Density metrics were computed as the proportion of echoes at different height percentiles relative to all echoes (p1, . . ., p99). In addition, mean (havg) and standard deviation (hstd) of the heights and the proportion of vegetation (veg) hits vs ground hits using thresholds of 0.5, 2, and 5 m were calculated. The following intensity metrics were also computed: percentiles of 10, 30, 50, 70 and 90 (i10 . . ., i90), the average intensity (iavg) and the standard deviations (istd). Some FUSION cloud metrics (McGaughey, 2012), such as skewness, kurtosis, AAD (average absolute deviation), L-moments (describes the characteristics of height distribution; mean, variance, skewness, kurtosis and their ratios) and elevation quadratic/cubic mean (mean of the height in power 2 and 3), were also calculated using first echo data. In addition, the proportion of first (fp pro), last (lp pro) and intermediate (ip pro) echoes from all echoes, the number of last echoes divided by the number of first echoes (lp fp), and the number of intermediate echoes divided by first echoes (ip fp) were calculated. Well-known ALS point cloud metrics were selected in this study due their strong correlation with canopy height and density. We expected the point cloud metrics to recognize different canopy characteristics especially between different LU/LC classes, site types and peat land types. ALS surface metrics (raster attributes) were calculated from ALSbased surface models (rasters), which were created from DTM, DSM and CHM models. GRASS GIS 6.4.3 (GRASS GIS, 2014) was used to generate slope, tangential (tcurv) and profile curvature (pcurv), wetness index (wet) and flowline density (fld) rasters, accumulation raster using single flow direction (asfd) and accumulation raster using multiple flow directions (amfd) calculated from DTM and DSM. To be able to create flowline density and accumulation rasters, surface conditions were taken into account from a 200 m buffer around the plot centre. Negative values of the wetness index were converted to zero, which is more logical minimum for this index to be used in the calculations and negative values of the accumulation rasters were changed to positive because this means that the runoff is coming outside of the calculation area. In addition, a binary tangential curvature raster (Ctcurv) and profile curvature raster (Cpcurv) were created from corresponding rasters by classifying convex values into one class and concave values into another class. QGIS 2.0.1 (QGIS, 2014) was used to create hillshade (hills) using the azimuth of 135◦ (northwest), 180◦ (north), 225◦ (northeast) and declination of 30◦ and ruggedness (rug, describes terrain heterogeneity) rasters from DTM and DSM. CHM was classified into two classes (C2 chm) using a height limit of 1 m, i.e. ground and vegetation pixels, and into seven classes (C7 chm) using height limits of 5, 10, 15, 20, 25 and 30 m. All rasters were created using a spatial resolution of 1 m. More detailed description of calculation of used surfaces/rasters can be found from manuals of used programmes.

16

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

Table 2 Point cloud, surface and spectral metrics used as predictor variables in classification models. Metrics Point cloud metrics f/l h1/2/10/20/. . ./90/95/100 f/l p1/5/10/20/. . ./90/95/99 f/l havg/hstd f/l veg05/2/5 f/l i10/30/50/70/90 f/l iavg/istd Elev CV/IQa Elev skewnessa Elev kurtosisa Elev AADa Elev MAD median/modea Elev L1/L2/L3/L4a canopy relief ratioa Elev L CV/skewness/kurtosisa Elev SQRTm2/CURTm3a fp/lp/ip pro lp/ip fp Surface metrics a/s dtm/dsm/chm a/s slope dtm/dsm a/s pcurv dtm/dsm a/s tcurv dtm/dsm a/s wet dtm/dsm a/s fld dtm/dsm a/s asfd dtm/dsm a/s amfd dtm/dsm a/s hills135/180/225 dtm/dsm a/s rugg dtm/dsm a/s Cpcurv dtm/dsm a/s Ctcurv dtm/dsm ASM/cont/cor/texvar/idm/savg/svar/sentro/entro/dvar/dent 10/30 Cpcurv dtm/dsm

Description Canopy height percentiles Proportional canopy densities of cumulative heights Average and standard deviation of echo heights Proportion of vegetation hits above ground, vegetation limit 0.5, 2 and 5 m Intensity percentiles of cumulative heights Average and standard deviation of intensity values Coefficient of variation and interquartile distance of heights Skewness (describes asymmetry of distribution) Kurtosis (describes peaky of distribution) Average absolute deviation of heights Median and mode of the absolute deviation of heights from overall median and mode L-moments (describes the characteristics of distribution) Canopy relief ratio L-moment coefficient variation, skewness and kurtosis Quadratic and cubic mean of heights (power 2 and 3) Proportion of echo types (first, last, intermediate) Ratio between echo types (last/first, intermediate/first)

ASM/cont/cor/texvar/idm/savg/svar/sentro/entro/dvar/dent C7 chm

Average and standard deviation value of DTM and DSM and CHM raster Average and standard deviation value of DTM and DSM raster As above but for profile curvature raster As above but for tangential curvature raster As above but for wetness index raster As above but for flowline density raster As above but for accumulation raster (single flow direction) As above but for accumulation raster (multiple flow direction) As above but for hillshade raster with different azimuth As above but for ruggedness raster As above but for classified profile curvature raster As above but for classified tangential curvature raster Haralick et al.’s texture metrics (1973) using requantification classes of 10 and 30 for the classified profile curvature raster calculated from DTM and DSM As above but for tangential curvature raster Average and standard deviation value of classified canopy height model (two classes) calculated from CHM As above but for seven classes Haralick’s texture metrics (1973) for classified canopy height model (two classes) calculated from CHM As above but for seven classes

Spectral metrics a/min b1/b2/b3/b4/b5/b6/b7 b1/b2,. . .,b6/b7

Average and minimum value of band 1–7 Ratio between average values of two different bands 1–7

ASM/cont/cor/texvar/idm/savg/svar/sentro/entro/dvar/dent 10/30 Ctcurv dtm/dsm a/s C2 chm a/s C7 chm ASM/cont/cor/texvar/idm/savg/svar/sentro/entro/dvar/dent C2 chm

Prefixes: f = first echo data, l = last echo data, s = standard deviation, a = average. The suffixes DTM, DSM, CHM show the surface model from which the metrics have been calculated. a Point cloud metrics calculated using FUSION software.

Following this, the average (a) and standard (s) deviations of the pixel values in all the rasters (surface models) were calculated inside the boundaries of the circular plots and segments. Furthermore, Haralick et al.’s texture metrics (1973) were also calculated from the classified curvature rasters and classified CHM rasters. Texture metrics were calculated by creating only one grey tone spatial-dependence matrix from the area of a plot or segment. An average of all directions (0, 45, 90 and 135◦ ) with a lag value of 1 m was used. Lag value of 1 m was selected after preliminary tests. It corresponds to the lag of one pixel. ALS-based surface metrics were expected to better describe the forest land attributes than straightforward point cloud metrics. This is because forest land attributes, such as main type, site type, peat land type and drainage status are strongly related to the properties of the ground/soil, which can be described by DTM based metrics such as slope and ruggedness. On the other hand, DSM and CHM based metrics were also tested, since they could also describe different properties of surface/canopy such as variation, cover and shadowness compared to point cloud metrics. They were expected to be especially useful in classification of LU/LC classes, site types and peat land types.

Spectral metrics were calculated from Landsat satellite images. The spectral metrics contained average (a) and minimum (min) values of each of seven bands and ratios between averages (e.g. b4b5 = average band 4/average band 5) inside the plots and segments. Based on earlier studies spectral metrics were expected to detect differences in tree species and vegetation cover and would, thus, be useful for classification of LU/LC classes, site type, peat land type. 2.5. Classification and variable selection Multinomial logistic regression analysis was used to classify plots and segments by means of the point cloud, surface and spectral metrics. Multinomial logistic regression is used to predict the probabilities of different possible outcomes of a categorically distributed dependent variable, given a set of independent variables (e.g. Greene, 2002). In the multinomial logit model, log-odds of each response are assumed to follow a linear model: Xij = log

pij piJ

= ˛j + xi ˇj

(1)

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

17

where i is a plot, j is a category and J is a number of categories, ˛j is a constant and ˇj is a vector of the regression coefficient, for j = 1,2, . . ., J − 1. Here constant is written explicitly, so we assume henceforth that the model matrix X does not include a column of ones. The multinomial logit model can also been written in terms of original probabilities pij by starting with Eq. (1) and adopting the convention that piJ = 0: pij =

J

exp{Xij }

k=1

exp{Xik }

(2)

for j = 1, . . ., J. Both j and k are the categories. The model ensures that all probabilities are in the interval [0,1] and that the probabilities sum to 1. The model coefficients were estimated with the method of maximum likelihood in R (www.r-project.org) using the mlogitfunction provided by Croissant (2014). Final classification was done by converting the highest probabilities into corresponding classes. Classification accuracies were tested by means of overall accuracy (OA, percent of all observations classified correctly), producers and user’s accuracy (PA and UA) and well-known kappa-value (kappa) (Landis and Koch, 1977). First, classification models were estimated based on the full data separately for circular plots and segments. We also tested the approach where training data were divided in two sub-data sets (stratified data), and the classification model was estimated separately for both sub-data sets. The metric f veg2 (proportion of vegetation hits vs all hits above 2 m height) was found to be the best performing metric for stratifying the data. Sub-data 1 consisted of the plots or segments where the value of f veg2 was 50% or more and sub-data 2 consisted of the rest of the plots or segments. By stratification we expected that ALS metrics from plots with a high proportion of either ground or vegetation hits would have a different relationship with the classes of interest. Variable selection for all models was done in two phases. In the first phase each metric was used as a single predictor variable in classification model. Variables were ordered based on the kappa-values and OAs of the models and circa 50 of the best performing metrics were selected for the second phase. When metrics were listed, some metrics obtained the same kappa-values and all those metrics were selected even if total number of metrics then exceeded 50. In addition, at least one best metric from each variable group (point cloud, surface, spectral) was added to the list of best performing metrics. This was done to increase the variability of metrics, which are not so correlated with each other. In the second phase all possible combinations of the 50 best performing metrics were tested. The maximum number of predictor variables in the classification models was fixed at three to avoid overfitting. The best performing metrics for models were selected based on the kappa-value of classifications and correlations among metrics. Final classification accuracies were calculated based on leave-one-cluster-out cross-validation, which means that the whole NFI cluster (max 9–12 plots) in which the target plot or segment belonged was excluded when the model was fitted. Workflow of classification approach is shown in Fig. 4. 3. Results

Fig. 4. Workflow of classification approach. Table 3 Overall accuracies (OA) and kappa-values (kappa) of the most accurate classification models using plots (P) and segments (S) and full and stratified data. Classification

Full

National LU/LC, P National LU/LC, S FAO LU/LC, P FAO LU/LC, S Main type, P Main type, S Drainage, P Drainage, S Site type, P Site type, S Peat land type, P Peat land type, S

Stratified

OA

Kappa

OA

Kappa

0.94 0.96 0.95 0.95 0.90 0.90 0.86 0.88 0.68 0.69 0.81 0.82

0.80 0.84 0.82 0.82 0.24 0.34 0.40 0.52 0.48 0.50 0.60 0.63

0.94 0.97 0.96 0.97 0.90 0.90 0.86 0.88 0.68 0.69 0.84 0.84

0.80 0.90 0.86 0.91 0.32 0.37 0.38 0.49 0.48 0.51 0.69 0.69

Table 4 Confusion matrix of the national LU/LC classification. Forestry land

Agr. land

Built areas

Total

UA

Forestry land Agr. land Built areas

673 6 5

4 104 1

6 3 26

683 113 32

99% 92% 81%

Total PA

684 98%

109 95%

35 74%

828

PA = producer’s accuracy, UA = user’s accuracy, Agr. = agricultural.

Table 5 Confusion matrix of the FAO LU/LC classification. Forest

Non-forest

Total

UA

Forest Non-forest

665 10

12 141

677 151

98% 93%

Total PA

675 99%

153 92%

828

3.1. Classification accuracies Classification of the national and FAO LU/LC classes resulted in the highest OAs and kappa-values in this study (OA = 0.97, kappa = 0.90 and OA = 0.97, kappa = 0.91) (Table 3). Forestry land and agricultural land were classified more accurately than class built and forest slightly better than non-forest (Tables 4 and 5). The main type and drainage status also obtained high OAs (0.90 and 0.88), but the kappa-values of these classifications were only 0.37

and 0.52. Almost all mineral soils and undrained areas were classified correctly, but only 30% of peat lands and half of drained areas (Tables 6 and 7). Classification of the site types resulted in an OA of 0.69 and kappa of 0.51. The medium site type was classified slightly better than rich and poor (Table 8). The most accurate classification model for peat land type gave an OA of 0.84 and kappa of 0.69 (Table 9). Spruce peat lands were classified the most accurately.

18

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

Table 6 Confusion matrix of the main type classification. Mineral soil

Peat land

Total

UA

Mineral soil Peat land

587 17

54 23

641 40

92% 58%

Total PA

604 97%

77 30%

681

Table 7 Confusion matrix of the drainage status classification. Drained Drained Undrained

57 60 117 49%

Total PA

Undrained

Total

UA

21 453

78 603

73% 90%

564 96%

681

Rich Medium Poor Total PA

89 4 37 130 69%

Medium

Poor

Total

UA

3 192 71

29 66 190

121 262 298

74% 73% 64%

266 72%

285 67%

681

Open

Total PA

3 1 0 4 75%

Metrics

National LU/LC, P National LU/LC, S FAO LU/LC, P FAO LU/LC, S Main type, P Main type, S Site type, P Site type, S Peat land type, P Peat land type, S Drainage, P Drainage, S

s hills225 dsm s hills225 dsm s rtcurv dsm s rtcurv dsm a slope dtm s rpcurv dtm b4b5 b4b5 l p5 f veg2 a rrug dtm s rpcurv dtm

s slope dsm s rtcurv dtm s rtcurv dtm ElevL2a l iavg a rugg dtm a b7 b5b6 ElevLskewnessa l veg05 s rpcurv dtm a rugg dtm

l istd l istd f istd l istd s slope dtm a amfd dtm l hstd ip pro a C2 chm s wet dsm a C2 chm a asfd dtm

both full and stratified data. Table 10 shows an example of used metrics in the most accurate classification models using full data. 4. Discussion 4.1. Evaluation of chosen metrics and combined use of ALS and satellite data

Table 9 Confusion matrix of the peat land type classification.

Open Spruce Pine

Classification

Prefixes: f = first echo data, l = last echo data, s = standard deviation, a = average. The suffixes DTM, DSM and CHM show the base surface where the metrics have been calculated. a Point cloud metrics calculated using FUSION software.

Table 8 Confusion matrix of the site type classification. Rich

Table 10 Predictor variables (point cloud, surface and spectral metrics) used in the best performing classification models using full data and plots (P) and segments (S).

Spruce

Pine

Total

UA

2 43 3

0 6 19

5 50 22

60% 86% 86%

48 90%

25 76%

77

In most classification schemes the use of segments gave slightly better kappa-values than the use of circular plots. In particular the classification of the main type and drainage status improved considerably. Only small differences were observed in classification accuracies between use of full data and stratified data. In particular, classification of FAO LU/LC class, peat land type and main type succeeded slightly better using stratified data compared to full data. Also use of segments in stratified data gave higher accuracies for most classification schemes compared to use of plots (Table 3). In general, the incorrectly classified plots and segments were classified into the closest or most similar class, e.g. site type poor into medium and LU/LC class built into forest land and peat land type pine into type spruce. Classes with a low number of observations were generally classified more often incorrectly than classes with a high number of observations. 3.2. Best performing metrics In the classification of national and FAO LU/LC class, main type and drainage status the surface metrics were found to be the most important predictor variables. Point cloud metrics, especially intensity metrics, were used together with surface metrics in the most accurate models. Best performing models for main type and drainage status favoured DTM-based surface metrics. Spectral metrics were found to be the most important metrics in the case of site type. Different point cloud metrics were used together with spectral metrics. Classification of peat land types was only case where point cloud metrics were found to be the most important metrics. Similar metrics were found to be important variables using

The reason for the successful use of ALS-based surface metrics calculated from DSM in classification of LU/LC classes might be that they detect better the horizontal properties such as variation, steepness and cover of the surface/canopy compared to other metrics. Point cloud metrics describe the vertical structure of canopy and spectral metrics the spectral value of target. Neither point cloud nor spectral metrics separated different LU/LC classes as well as surface metrics. In the classification of main type and drainage status surface metrics calculated from DTM were the most important, which shows that properties of ground/soil, such as variation in terrain and water conditions, are more important when separating mineral soil and peat land and drained and undrained classes than properties of canopy (point cloud metrics) or spectral value (spectral metrics). This study also showed that site fertility type can be most accurately detected based on spectral values (spectral metrics), because they recognize the proportion of deciduous trees and other green vegetation in different site types. ALS-based metrics cannot detect differences in ground/soil and canopy structure (surface and point cloud metrics) accurately enough to separate site fertility types compared to spectral metrics, even though some point cloud metrics were used together with spectral metrics. The reason for the successful use of point cloud metrics (and some DSM-based surface metrics) in classification of peat land types might be that they are related to the differences in canopy densities between open peat land, pine dominated peat land and spruce dominated peat land. DTM-based surface metrics did not perform correspondingly because general properties of ground/soil are quite similar in different types of peat lands and spectral metrics cannot detect canopy structure in detail. This study showed that it is worth to use remote sensing data from multiple data sources because different sensors detect different characteristics. Forest site type was most accurately detected using spectral metrics from satellite images. Correspondingly, LU/LC classes, main type, peat land type and drainage status were detected using point cloud and surface metrics from ALS data. It is also profitable to use metrics from the different data sources together in the same model, as it was done in the case of site type when point cloud metrics were used together with spectral metrics.

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

4.2. National and FAO LU/LC class Classification of LU/LC class according to both national and FAO classification was highly accurate and ALS-based surface metrics, especially calculated from DSM, played a major role as predictor variables. For example, shadow, slope and curvature conditions of the canopy in the forestry land or forest class had much more variation compared to the agricultural land and non-forest categories. Traditionally these surfaces are calculated from the DTM. Some ALS intensity metrics, which produced a higher deviation of intensity values in forestry land and forest, were also found to be good predictors. Significance of surface and intensity metrics is reasonable, because structure and reflectance of canopy surface between forestry land, agricultural land and built-up areas or forest and non-forest differ. Misclassification of built and non-forest classes can be explained by the heterogeneity within these classes. Those classes included, e.g., built-up areas, power lines, wooded lots, roads and other wooded land which all have quite different surface conditions. Some of those have very similar characteristics to that of forestry land, which is one reason for misclassifications. Our high classification accuracies of forest land and forest are in line with earlier land use classification studies where ALS-based surfaces have been utilized. For example, Charaniya et al. (2004) obtained an OA of 66–84% using supervised parametric classification. They discriminated roads, grass, roofs and trees and used height, height variation, height difference between first and last returns and intensity rasters in their classification. Brennan and Webster (2006) used DSM, DEM, height, intensity and multiple return rasters in object-oriented classification and separated data into ten classes with an OA and kappa of 94% and 0.93. Antonarakis et al. (2008) also used object-based classification for nine classes with an OA of 95%. They used vegetation height, percentage canopy hits, intensity, skewness and kurtosis rasters and classified five forest types with accuracies between 66 and 98%. Earlier studies are not straightforwardly comparable with this study due to quite different class definitions. 4.3. Main type and drainage status This was the first time when mineral soils and peat lands and undrained and drained forest land was separated from each other using ALS data. Although the classification accuracies of main type and drainage status were only low and moderate, we found that ALS-based surface metrics, especially calculated from DTM, were important predictors. We obtained, for example, a higher average or deviation in slope, curvature and ruggedness rasters in the mineral soil and undrained classes compared to the peat land and drained classes. This might be explained by the flat and smooth terrain surface of peat lands and drained areas compared with mineral soils and undrained areas. Some hydrological surface metrics were also found to have rather good predictive power, which was expected because of the higher water accumulation of peat lands (flat) and drained areas (ditches) compared with mineral soil and undrained areas. Misclassification rates of peat lands and drained areas were high. This means that classes peat land and drained areas cannot be detected accurately enough using ALS and satellite data. One reason for the misclassification of peat lands may be that large proportions of the peat lands in the study area have already changed very close to mineral soils because of the successful draining. Some characteristics of drained areas might not be visible because of overgrown ditches or a sparse ditch network. Fig. 5 shows an example of two drained plots/segments. In the left image ditches cannot be seen from DTM and the plot/segment is misclassified, but in the right image they are clearly visible and, thus, correctly classified. Differences (heterogeneity) among peat lands and among drained areas might also affect the difficulty of separating between mineral soils

19

and peat lands and undrained and drained areas. In addition, the number of observations in classes peat land and drained areas was relatively small compared to classes mineral soils and undrained, which also might have affected to the results. We could not find any ALS-based studies where mineral soils and peat lands were separated from each other, but Maxa and Bolstad (2009) discriminated upland, water and eight wetland classes using satellite images and ALS data with a kappa value of 0.7. ALS-based surfaces used in the classification were: elevation, slope and terrain shape. Dirksen (2013) discriminated swamp forest and non-paludified forests in old-growth forests using ALS point cloud metrics. According to their results vegetation height is lower in swamp forests than in non-paludified forests. However, the accuracy was only 54–62% due to a border effect of small swamps and differences in structure among and within individual swamps. 4.4. Site type This study showed that combination of spectral and point cloud metrics can be used to predict forest site types. Spectral metrics calculated from satellite images played a more important role in the classification of site types than ALS-based metrics. Especially metric b4b5 was very important. Metric b4b5 produced higher values the richer the site type was. This indicates that a combination of band 4 (near infrared) and band 5 (short wave infrared) recognizes the larger amount of deciduous trees and green vegetation in more fertile site types, which is to be expected. From point cloud metrics especially the proportion of intermediate echoes was often used, but also some height-related metrics appeared in models. Those metrics produced higher values the richer the site type was. Vehmas et al. (2011) also found higher amount of intermediate echoes in more fertile sites, which can be explained by a larger amount of understory vegetation, resulting in a larger number of intermediate echoes. The reason for the behaviour of the height metrics is probably related to the fact that trees grow taller in fertile site types. So far forest site types have been classified using only one remotely sensed data-source, such as satellite data or ALS data. The significance of spectral metrics calculated from satellite and aerial images for the classification of forest types has been shown in many studies (e.g. Tomppo, 1992; Tomppo et al., 2009). Tomppo (1992) classified pixels into four site type classes for forest income taxation purposes using Landsat TM satellite images. Their pixel level OA varied between 50% and 70%, and the richest and poorest site types were classified most accurately. Furthermore, in their study the ratios between bands worked better than the original bands and similarly the best performing bands were band 4 (near infrared) followed by bands 7 and 1. So far only a couple of site-type classification studies have utilized ALS data. Vehmas et al. (2011) classified five forest site types with a kappa value of 0.47 using ALS height and intensity metrics as predictor variables. Even though their field data consisted of 247 forest stands in mineral soil of a mature age, their classification accuracies were similar to our results. Holopainen et al. (2010) first used dominant height predicted from ALS data and stand age from stand register to determine the site index and then converted site indexes into site types. They got a kappa value of 0.60 for five site types. Vehmas et al. (2009) identified herb-rich forests from less fertile site types using ALS point cloud variables and a logistic and k-NN classifier with a kappa value of 0.64. They used the same data as Vehmas et al. (2011) and the classification accuracy was 65% for herb-rich stands and 96% for other stands. 4.5. Peat land type It was only in the classification of peat land types that the point cloud metrics, especially density metrics, were found to be

20

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

Fig. 5. Example of two drained plots/segments. In the left image e plot/segment is incorrectly classified into class undrained and ditches cannot be seen from DTM. In the right image ditches are clearly visible and plot/segment is correctly classified.

most significant predictor variables. For example, spruce peat lands obtained a higher proportion of vegetation hits than pine and open peat lands. This can be expected because canopy cover in spruce peat lands is usually higher than in pine peat lands and open peat lands, which have very low canopy cover. Some ALS-based surface metrics calculated from CHM (such as average value of classified canopy height model (two classes)) and DSM (such as standard deviation of wetness index raster) were also selected. They indicated the highest canopy cover as being spruce peat land, then pine peat land and the lowest to open peat land and a higher variation in surface level in spruce and pine peat land compared with open peat land, as expected. Korpela et al. (2009) tested ALS-based features in the classification of mire habitats. They found that surface characteristics of mire types could be delineated with high accuracy, but qualifying differences in ground flora failed. Their 21 mire types were classified with a kappa of 0.25–0.62 and corresponding values for the main mire type (treeless, composite and forested) were 0.32–0.67 and for dominant tree species (spruce, pine and spruce–pine) 0.66–0.81. Their classification accuracies of the main mire types and dominant tree species are in line with our results, but are not fully comparable, because our classification consisted of pine dominated, spruce dominated and open peat land classes. They also found point cloud metrics (point proportions and height distribution metrics) to be the most significant predictor variables in classification models and also intensity metrics in the separation of tree species, even though few variables calculated from DTM were tested. 4.6. Reference plot size and stratification The use of segments gave slightly more accurate results in almost all classification schemes compared with the use of circular plots. In particular, the classification of main type and drainage status, which utilized mainly DTM-based metrics, benefited significantly from larger segments. The reason for this may be that larger segments represent the characteristics of forest land and site type of the stand better than smaller plots. Larger segments reduce the random variation inside segments and, e.g. characteristics of sparse ditch network are then more representative (see Fig. 4, right image). Hence, forest land can be classified based on terrain characteristics, such as slope, curvature, ruggedness and water accumulation, if metrics have been calculated from an area large enough. Classification of site type benefited less from the use of segments, which is reasonable, because mainly spectral metrics were used in best performing classification models of site types. Use of the stratified data (stratification based on proportion of vegetation hits) did not

improve the classification accuracies significantly. The reason for this could be that many of the classifications did not use forest canopy structure metrics (e.g. point density metrics and surface metrics calculated from DSM and CHM) in classification models. Those classifications which did (e.g. FAO LU/LC class and peat land type) obtained slightly better results. 4.7. Applying forest land attributes in practical inventory Each ALS inventory usually applies new remote sensing data and local field training data always need to be collected to ensure the success of the inventory. Hence, in different inventory areas different metrics might be useful, e.g. if tree species proportions or properties of forest and forest land differ. Training data should cover all the variation in an inventory area, such as all considered site types in different developing stages and dominant tree species. In addition, field data and remote sensing data should be acquired during the same growing period, even though forest land attributes do not change as fast as growing stock attributes. If major changes, such as clear cuts have happened, between the acquisition of field data and remote sensing data, those samples should be removed from the training data. At the moment, ALS-based forest inventory applies point cloud metrics together with spectral metrics from aerial images to predict species specific growing stock attributes for large areas (Maltamo & Packalen 2014). Prediction of site type and forest land categories could be incorporated by including registration of forest land attributes during field work, adding satellite images and calculating ALS-based surface metrics from ALS point cloud. Combined use of remote sensing data and information from existing land and stand register data might also be applied in practical inventories in the future (Pippuri et al., 2013). In practical inventory, when predicting the forest land categories and site types, similar segments from CHM, such as the ones we used in this study, could be delineated using automatic segmentation. However, because our results were only slightly better using segments in most of the classification schemes, the gridbased approach could also be applicable (Fig. 6). Hence, in this study the plot level modelling unit was large enough to describe forest land characteristics in most cases. At the moment in ALSbased inventories prediction is done in grid cells or homogeneous micro-stands. In both cases forest land and site type information can be generalized to the stand level or, e.g., used for stand delineation. Unfortunately, both grid cells and micro-stands can include more than one forest land or site type classes; micro-stands are usually formed based on tree dimensions, not site type.

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

21

Novel applications for these metrics and combined use of data sources might also be found from different areas of forest sciences. Acknowledgements This work was funded by the doctoral student position and strategic funding of University of Eastern Finland, which are acknowledged. References

Fig. 6. Example of prediction of site fertility types (Cajander, 1926) using a grid approach. Site type rich = Oxalis–Maianthemum and Oxalis–Myrtillus types, medium = Myrtillus type, poor = Vaccinum, Calluna and Cladina types and rock and sand sites.

NFI field measurement as training data was successfully utilized in this study. Using field plots from the NFI, would increase the costeffectiveness of ALS-based inventory, since separate field data do not need to be collected. In addition, the accuracy of the inventory will improve when classification of forest land attributes are done by an NFI field expert. Earlier, e.g. Tuominen et al. (2014) provided promising results for growing stock attributes using NFI plot data in ALS-based forest inventories. Inventory method presented in this study is an accurate and cost-effective way to predict forest land attributes for large areas. If more accurate results are needed, it will require field inventory. The comparison of the results of this study to earlier works is difficult due to the different classification cases, methods and study areas.

5. Conclusions This study showed that forest land attributes can be classified using ALS data, satellite images and field plots from National Forest Inventory of Finland. Prediction of site type and forest land category could be incorporated into a stand level forest management inventory system in Finland, as presented in this study. Both ALS data and satellite data were applied in the classifications of forest land attributes, which supports the benefit of using multiple sources of remotely sensed data in forest inventories. Rarely used ALS-based surface metrics were found to be useful predictor variables in classification of land use/cover, main type and drainage status. They take into account the horizontal properties of surface/canopy and ground/soil in more detail compared to traditional ALS point cloud metrics and spectral metrics. Spectral metrics were found to be the most useful metrics in the classification of site type whereas ALS point cloud metrics were the most useful metrics in the classification of peat land type. Based on the findings of this study ALS-based surface metrics and combined use of multiple data sources should be utilized in practical forestry to develop more costeffective forest inventory and management planning applications.

Antonarakis, A.S., Richards, K.S., Brasington, J., 2008. Object-based land cover classification using airborne LiDAR. Remote Sens. Environ. 112 (6), 2988–2998. Axelsson, P., 2000. DEM generation from laser scanning data using adaptive TIN models. Int. Arch. Photogramm. Remote Sens. 33 (Part B4), 110–117. Brennan, R., Webster, T.L., 2006. Object-oriented land cover classification of Lidarderived surfaces. Can. J. Remote Sens. 32 (2), 162–172. Cajander, A.K., 1926. The theory of forests types. Acta For. Fenn. 29 (3), 1–108. Charaniya, A.P., Manduchi, R., Lodha, S.K., 2004. Supervised parametric classification of aerial LiDAR data. In: Proc. 2004 IEEE Computer Society Conference on Computer Vision and Pattern Recognition Workshops (CVPRW’04), vol. 3, Washington, DC, USA, 27 June–2 July 2004, p. 30. Croissant, Y., 2014. Package ‘mlogit’, http://cran.r-project.org/web/packages/mlogit/ mlogit.pdf Dirksen, J.W., 2013. Modelling Presence of Swamp Forest and Forest Dwelling Birds in a Boreal Forest Reserve Using Airborne Laser Scanning. Norwegian University of Life Sciences, Department of Ecology and Natural Resource Management, Ås. FAO, 2006. Global Forest Resources Assessment 2005. Progress Towards Sustainable Forest Management. Food and Agricultural Organization of the United Nations (FAO), Rome, FAO Forestry paper no. 147. Gatziolis, D., 2007. Lidar-derived site index in the U.S. Pacific Northwest – challenges and opportunities. The international archives of the photogrammetry. Remote Sens. Spat. Inf. Sci. 36 (Part 3), 136–143. GRASS GIS, 2014. www.grass.osgeo.org (accessed 28.02.14). Greene, W.H., 2002. Econometric Analysis, fifth ed. Prentice Hall, Upper Saddle River, NJ. Haralick, R.M., Shanmugam, K., Dinstein, I., 1973. Textural features for image classification. IEEE Trans. Syst. Man Cybern. 3 (6), 610–621. Hollaus, M., Wagner, W., Maier, B., Schadauer, K., 2007. Airborne laser scanning of forest stem volume in a mountainous environment. Sensors 7, 1559–1577. Holmgren, J., 2004. Prediction of tree height, basal area and stem volume in forest stands using airborne laser scanning. Scand. J. For. Res. 19 (6), 543–553. Holopainen, M., Vastaranta, M., Haapanen, R., Yu, X.W., Hyyppä, J., Kaartinen, H., Viitala, R., Hyyppä, H., 2010. Site-type estimation using airborne laser scanning and stand register data. Photogramm. J. Finl. 22 (1), 16–32. Hyyppä, J., Inkinen, M., 1999. Detecting and estimating attributes for single trees using laser scanner. Photogramm. J. Finl. 16 (2), 27–42. Koch, B., Straub, C., Dees, M., Wang, Y., Weinacker, H., 2009. Airborne laser data for stand delineation and information extraction. Int. J. Remote Sens. 30 (4), 935–963. Korhonen, R., Korpela, L., Sarkkola, S., 2008. Finland – Fenland. Research and Sustainable Utilisation of Mires and Peat. Finnish Peatland Society and Maahenki Ltd., Helsinki. Korhonen, L., Korpela, I., Heiskanen, J., Maltamo, M., 2011. Airborne discrete-return LIDAR data in the estimation of vertical canopy cover, angular canopy closure and leaf area index. Remote Sens. Environ. 115 (4), 1065–1080. Korpela, I., Koskinen, M., Vasander, H., Holopainen, M., Minkkinen, K., 2009. Airborne small-footprint discrete-return LiDAR data in the assessment of boreal mire surface patterns, vegetation, and habitats. For. Ecol. Manag. 258 (7), 1549–1566. Laine, J., Vasander, H., 2008. Suotyypit ja niiden tunnistaminen, second ed. Metsäkustannus Oy, Helsinki (in Finnish). Landis, J.R., Koch, G.G., 1977. The measurement of observer agreement for categorical data. Biometrics 33 (1), 159–174. Maltamo, M., Packalen, P., 2014. Species-specific management inventory in Finland. In: Maltamo, M., Næsset, E., Vauhkonen, J. (Eds.), Forestry Applications of Airborne Laser Scanning: Concepts and Case Studies. Springer, Dordrecht, pp. 241–252. Maltamo, M., Packalén, P., Suvanto, A., Korhonen, K.T., Mehtätalo, L., Hyvönen, P., 2009. Combining ALS and NFI training data for forest management planning: a case study in Kuortane, Western Finland. Eur. J. For. Res. 128 (3), 305–317. Maxa, M., Bolstad, P., 2009. Mapping northern wetlands with high resolution satellite images and LiDAR. Wetlands 29 (1), 248–260. McGaughey, R.J., 2012. FUSION/LDV: Software for LIDAR Data Analysis and Visualization, May 2012 – FUSION Version 3.10. United States Department of Agriculture, Washington, DC. Næsset, E., 2002. Predicting forest stand characteristics with airborne scanning laser using a practical two-stage procedure and field data. Remote Sens. Environ. 80 (1), 88–99. Næsset, E., 2004. Practical large-scale forest stand inventory using small-footprint airborne scanning laser. Scand. J. For. Res. 19 (2), 64–179. Packalén, P., Maltamo, M., 2007. The k-MSN method for the prediction of speciesspecific stand attributes using airborne laser scanning and aerial photographs. Remote Sens. Environ. 109 (33), 328–341.

22

I. Pippuri et al. / International Journal of Applied Earth Observation and Geoinformation 44 (2016) 11–22

Packalén, P., Mehtätalo, L., Maltamo, M., 2011. Als-based estimation of plot volume and site index in a eucalyptus plantation with a nonlinear mixed-effect model that accounts for the clone effect. Ann. For. Sci. 68 (6), 1085–1092. Pippuri, I., Kallio, E., Maltamo, M., Peltola, H., Packalén, P., 2012. Exploring horizontal area-based metrics to discriminate the spatial pattern of trees and need for first thinning using airborne laser scanning. Forestry 85 (2), 305–314. Pippuri, I., Maltamo, M., Packalen, P., Mäkitalo, J., 2013. Predicting species-specific basal areas in urban forests using airborne laser scanning and existing stand register data. Eur. J. For. Res. 132 (5–6), 999–1012. QGIS, 2014. http://www.qgis.org (accessed 28.02.14). Racine, E.B., Coops, N.C., St-Onge, B., Bégin, J., 2014. Estimating forest age from LiDARderived predictors and nearest neighbor imputation. For. Sci. 60 (1), 128–136. Sader, S.A., Ahl, D., Liou, W.-S., 1995. Accuracy of Landsat-TM and GIS rule-based methods for forest wetland classification in Maine. Remote Sens. Environ. 53 (3), 133–144. Tokola, T., Löfman, S., Erkkilä, A., 1999. Relative calibration of multitemporal Landsat data for forest area change detection. Remote Sens. Environ. 68 (1), 1–11. Tomppo, E., 1992. Satellite image-aided forest site fertility estimation for forest income taxation purposes. Acta For. Fenn. 229, 1–70. Tomppo, E.O., Gabliano, C., Natale, F., Katila, M., McRoberts, R., 2009. Predicting categorial forest variables using an improved k-neasrest neighbour estimator and Landsat imagery. Remote Sens. Environ. 113, 500–517.

Townsend, P.A., Walsh, S.J., 2001. Remote sensing of forested wetlands: application of multitemporal and multispectral satellite imagery to determine plant community composition and structure in southeastern USA. Plant Ecol. 157 (2), 129–149. Tuominen, S., Pitkänen, J., Balazs, A., Korhonen, K.T., Hyvönen, P., Muinonen, E., 2014. NFI plots as complementary reference data in forest inventory based on airborne laser scanning and aerial photography in Finland. Silva Fenn. 48 (2), http://dx. doi.org/10.14214/sf.983 Van Aardt, J.A.N., Wynne, R.H., Scrivani, J.A., 2008. Lidar-based mapping of forest volume and biomass by taxonomic group using structural homogeneous segments. Photogramm. Eng. Remote Sens. 74 (8), 1033–1044. Vastaranta, M., Korpela, I., Uotila, A., Hovi, A., Holopainen, M., 2012. Mapping of snow-damaged trees based on bitemporal LiDAR data. Eur. J. For. Res. 131 (4), 1217–1228. Vehmas, M., Eerikäinen, K., Peuhkurinen, J., Packalén, P., Maltamo, M., 2009. Identification of boreal forest stands with high herbaceous plant diversity using airborne laser scanning. For. Ecol. Manag. 257 (1), 46–53. Vehmas, M., Eerikäinen, K., Peuhkurinen, J., Packalén, P., Maltamo, M., 2011. Airborne laser scanning for the site type identification of mature boreal forest stands. Remote Sens. 3 (1), 100–116. Zhu, Z., Woodcock, C.E., 2012. Object-based cloud and cloud shadow detection in Landsat imagery. Remote Sens. Environ. 118, 83–94.