Abstract Multispectral LiDAR (light detection and ranging) data have been initially used for land cover classification. However, there are still high classification uncertainties, especially in urban areas, where objects are often mixed and confounded. This study investigated the efficiency of combining advanced statistical methods and LiDAR metrics derived from multispectral LiDAR data for improving land cover classification accuracy in urban areas. The study area is located in Oshawa, Ontario, Canada, on the Lake Ontario shoreline. Multispectral Optech Titan LiDAR data over the study area were acquired on 3 September 2014 in a single strip of 3 km2. Using the channels at 1,550 nm (C1), 1,064 nm (C2) and 532 nm (C3), LiDAR intensity data, normalized digital surface model (nDSM), pseudo normalized difference vegetation index (PseudoNDVI), morphological profiles (MP), and a novel hierarchical morphological profiles (HMP) were derived and used as features for the classification. A support vector machine classifier with a radial basis function (RBF) kernel was applied in the classification stage, where the optimal parameters for the classifier were selected by a grid search procedure. The combination of intensity, pseudoNDVI, nDSM and HMP resulted in the best land cover classification, with an overall accuracy of 93.28%.

Citation: Huo L-Z, Silva CA, Klauberg C, Mohan M, Zhao L-J, Tang P, et al. (2018) Supervised spatial classification of multispectral LiDAR data in urban areas. PLoS ONE 13(10): e0206185. https://doi.org/10.1371/journal.pone.0206185 Editor: Tayyab Ikram Shah, Western University, CANADA Received: June 15, 2018; Accepted: October 7, 2018; Published: October 24, 2018 This is an open access article, free of all copyright, and may be freely reproduced, distributed, transmitted, modified, built upon, or otherwise used by anyone for any lawful purpose. The work is made available under the Creative Commons CC0 public domain dedication. Data Availability: The data are owned by the Teledyne Optech (http://www.teledyneoptech.com) company. The authors confirm that they did not have any special access privileges to these data and other researchers can request the data in the same manner by contacting Paul LaRocque at Paul.Larocque@teledyneoptech.com. Funding: This work is jointly supported by the National Natural Science Foundation of China (Grant No. 41401421, 41701397) and the Youth Foundation of President of Institute of Remote Sensing and Digital Earth, Chinese Academy of Sciences (Y6SJ2400CX). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. Competing interests: The authors have declared that no competing interests exist.

Introduction Urban land cover mapping is important for urban land management and planning [1, 2]. Remote sensing technology, with a fast revisiting period and a large coverage, provides the primary source data for a better understanding of urban land cover [3]. To fully utilize the data, land cover maps are usually derived based on visual interpretation or automatic classification methods [4]. Regional and global urban areas are mainly monitored by mapping impervious surface areas, which have anthropogenic features through which water cannot infiltrate into the soil (e.g., roads, parking lots, and rooftops) [5]. Impervious surface areas have been monitored by various satellite sensors, such as the Operational Linescan System (OLS) from the Defense Meteorological Satellite Program (DMSP) that images nighttime light [6], 250–1000 m Moderate Resolution Imaging Spectroradiometer (MODIS) imagery [3,7], and 30 m Landsat imagery [8,9]. Even though the coarse to moderate spatial resolution of these image types describes the spatial extent of urban areas, they often are not able to resolve the extent and spatial arrangement of subtypes within urban areas (e.g., buildings, roads, and trees). High spatial resolution imagery (e.g., IKONOS and QuickBird) provides new opportunities to map detailed urban areas at a fine scale by providing more detailed observations of the Earth [10, 11]. Nonetheless, urban land cover classification is still a challenging task due to the spectral heterogeneity and structural diversity of the complex geospatial objects present [12]. Significant efforts have been made to advance image classification efficiency by focusing on object based image analysis [13–15] extracting spatial features [16, 17] and even scene based image analysis methods [18]. However, the effects of shadowing and relief displacement still pose considerable challenges [19]. Hence, active sensors such as airborne LiDAR (i.e., Light Detection and Ranging) data have been investigated for land cover classification in urban areas in the last decade [20]. Airborne LiDAR provides the 3-D coordinates of the survey area by collecting a cloud of laser range measurements [20]. Based on 3-D point clouds, LiDAR can be further interpolated into raster layers of surfaces, such as digital surface models (DSM), and intensity data [21]. To fully utilize the LiDAR data for land cover classification, several groups of methods have been proposed [20, 22–26]. The first category of methods directly uses LiDAR point clouds, e.g., LiDAR point clouds can be directly used for urban feature extraction [23]; or semantic point cloud interpretation can be performed based on optimal neighborhood selection [24]. The second category of methods mainly relies on LiDAR points derived products, intensity image [25] or both intensity and DSM image (i.e., height information) [26]. The study [26] done by Zhou suggests that LiDAR data alone (by combining intensity and DSM data) could potentially be useful to accurately map urban land cover. Another popular category of methods adopts the strategy to combine LiDAR data with optical images (e.g., multispectral or hyperspectral images) for detailed urban land cover mapping [22, 25, 27, 28]; DSM data compliment multi-spectral information from passive optical imagery to identify different land cover types. A survey of LiDAR data for urban land cover mapping can be found in [20]. However, most of the previous studies exploiting LiDAR for urban land cover mapping used single-band LiDAR due to limited data availability. Since backscattered energy from LiDAR depends on target materials, target roughness, and laser wavelength [29], single-band LiDAR has a restricted ability to discriminate land cover types. Hence, multispectral LiDAR holds greater promise to map urban land cover classes. New multispectral LiDAR data sensors (e.g., Multispectral Optech Titan LiDAR), which measures backscattered energies at different wavelengths, provide new opportunities to classify urban land cover effectively [30]. Since the release of the first commercial airborne multispectral LiDAR system, several studies have been tested to assess capabilities to produce more accurate land cover maps [29–32, 33–37]. For instance, Teo et al. demonstrated that multi-wavelength LiDAR can provide higher accuracy than single-wavelength LiDAR for land cover classification [29]. Bakula et al. [35] applied a maximum likelihood classifier in a six-class land cover classification and achieved an overall accuracy of 90%. A maximum likelihood classifier was also used to classify intensity and height images, and the authors found that height information is important for urban land cover mapping [34]. Fernandez-Diaz et al. [36] assessed capabilities of the Titan multispectral LiDAR for land cover classification, bathymetric mapping and canopy characterization. Zou et al. [37] adopted the object-based method and found that pseudo normalized difference vegetation index (pseudoNDVI) calculated from multispectral LiDAR may improve vegetation identification. In another study [33], an object-based analysis was also performed on multispectral airborne LiDAR data for land cover classification and map updating in Finland. Motivated by previous studies [29, 33–37], this research is focused on classification of the multispectral LiDAR intensity rasters, sinceraster format is a more convenient than point clouds for land cover mapping. Previous studies have shown that multispectral LiDAR performs better than single band LiDAR in land cover classification [29], and found pseudoNDVI [37,38] and DSM [34,35] to be useful for improving classification accuracy of multispectral LiDAR data. in the novelty of this studyis that we investigated the role of spatial information in improving urban land cover classification results employing multispectral LiDAR intensity data (i.e., spectral-spatial classification [39]). Spectral-spatial classification aims to improve classification accuracy by combining spatial information contained in the multispectral or hyperspectral data [39], to produce more accurate classification maps [40]. It is widely studied with a rich set of algorithm developments [39,41–43], especially in hyperspectral image analysis. Among these methods, morphological profiles (MP) perform well due to their ability to capture spatial relationships among different land cover types [41,42]. In mathematical morphology [44], two fundamental operators are erosion and dilation, which are applied to an image with a set of known shapes, called the structuring elements (SE). Different combinations of erosion and dilation constitute opening and closing operations, which are the building blocks of MP. Functionally, the opening operation can remove objects smaller than the structuring elements, while the closing operation can fill small holes and connect adjacent objects. The traditional MP operates on the characteristic image (mainly the first or first several principal components) of multispectral or hyperspectral images. Basically, it is performed assuming the image data are in the same vertical level; i.e., the morphological operation is performed only considering its spatial neighborhood pixels while lacking the capability to consider whether its neighborhood pixels obviously belong to another land cover class (e.g., the opening and closing operation for a tree pixel bordering grass pixels will impose effects on those grass pixels). However, the LiDAR-derived nDSM provides the vertical context of the image, thus offering opportunities to introduce vertical context while extracting the MP feature. To take full advantage of the nDSM data, a novel hierarchical morphological profiles (HMP) feature is proposed. Hence, the specific goals of this study are as follows: (1) to assess the usefulness of the MP for improving multispectral LiDAR data classification; and (2) to further study effectiveness of the proposed HMP.

Discussion Importance of spatial features for multispectral LiDAR classification In this rapidly changing world, timely and accurate classification of urban areas is crucial for urban planning and sustainable management [1, 55]. In this study, we are making use of cutting-edge multispectral LiDAR sensor (Optech Titan) technology and auxiliary information derived from the multispectral LiDAR for land cover classification; by combining multispectral LiDAR data and morphological profiles we classify terrains in urban areas. The classification results from our study demonstrate the capability of multispectral LiDAR data in recording a diversity of spectral signals from land cover objects and thereby underscores its potential for effective urban area classification. Herein, we classified the multispectral LiDAR intensity data into four classes, namely buildings, trees, roads, and grass, and obtained highest overall classification accuracies when a combination of intensity data, PseudoNDVI, nDSM, and HMP was used. Previous studies have shown that multispectral LiDAR data provides more benefits for urban land cover classification [29]; to improve classification accuracy, PseudoNDVI and DSM or nDSM features have been utilized [34, 37, 56–58]. Consistent with previous studies, this study further demonstrated that PseudoNDVI and nDSM can improve classification accuracy of a multispectral LiDAR intensity image. In addition, this study further demonstrated that inclusion of spatial features (MP and the proposed HMP) is markedly useful for deriving high accuracy classification results compared to methods employing only multispectral LiDAR intensity data; specifically, the road and building classes benefited the most in our case. In general, it is hard to classify the aforementioned classes due to multiple reasons; roads in the urban areas are filled with vehicles, thus resulting in highly reflective spots in the image (see Fig 6A); buildings are usually covered with different colored roofs (see Fig 5A and Fig 6A). Although the nDSM feature brings greatly improved classification accuracy for these two classes due to their distinct height distributions, it also introduced more classification errors between trees and grass (see the confusion matrix S2 Table), thus lowering the classification accuracies for these two classes. MP improved the classification accuracy of roads and buildings by 11.03% and 7.71%, respectively. However, many classification errors between these two classes resulted by including the MP feature (e.g., producer accuracy for the building class is only up to 50.13%), because it lacks the ability to operate in the vertical dimension. The proposed HMP considers the height information by extracting the MP feature over different vertical layers. In this way, the classification accuracy improved for these two classes by up to 98.23% and 93.01%. This demonstrates that the MP feature is useful for classification of multispectral LiDAR intensity data by capturing the spatial characteristics of different land cover types. By incorporating the vertical information (usually provided by the nDSM data from the LiDAR data), the proposed HMP feature provided better discrimination than the MP feature. Comparisons with previous studies The support vector machines classifier was adopted as the base classifier for this study due to its robustness and high accuracy for remote sensing image classification [54, 59, 60]. Stacking vectors of HMP with the multispectral LiDAR intensity image was found to increase classification accuracy significantly, resulting in an overall accuracy of 93.28% for the study area. Alternative methods that have been used for classification related to multispectral LiDAR data are object-based analysis, data clustering methods, and maximum likelihood classifiers. In Zou et al. [37], an Object Based Image Analysis (OBIA) approach for 3D land cover classification using multispectral LiDAR point clouds was presented; an overall accuracy of over 90% was obtained while classifying the landcover types into 9 categories. However, some misclassified objects–such as roads classified as lawn and bare soil–resulted due to similar spectral properties and at times because of lack of effective spectral and spatial features for distinguishing class boundaries. Matikainen et al. [33] evaluated the use of different spectral indices derived from multispectral LiDAR data for land cover classification and map updating; the classes considered for the study included building, tree, asphalt, gravel, rocky area and low vegetation, and obtained an accuracy of 96% compared with validation points. Multispectral Lidar data combined with old map vectors proved to enhance automated change detection of buildings, and assisted in removal of shadows on intensity images produced from the data; the passive aerial imaging commonly used in mapping suffered from external illumination conditions and often resulted in excessive shadows on intensity images. Bakula et al. [35] fused multi-wavelength laser intensity imagery, elevation data and textural information, and applied spectral (using maximum likelihood rule) and spectral-textural classification approaches for distinguishing 6 classes; they got an overall accuracy of 90%. In that study, applications of additional morphological classification and granulometric transformations were found to highly enhance the accuracy of the separation of building and road classes, as they eliminated several pixels that were initially confused. They also noticed that interpolating the intensity raster was not very helpful for improving classification results; even though using intensity rasters with both first and last returns slightly benefited the study. In essence, our strategy of using intensity features for classification resulted in accuracies similar to related studies, and considering the boost it gave to the overall classification (18.62%) obtained through the combination of intensity image with PseudoNDVI, nDSM, and HMP features, this method may be efficient for future multispectral LiDAR endeavors; e.g., plant species classification [61–63], urban change detection [64,65], flood inundation mapping [66] and even carbon sequestration modeling [67]. Implications and future directions Although multispectral LiDAR data can be considered a pivotal tool for bolstering subsequent urban planning and mapping operations [36, 38, 68], data processing should be undertaken with caution. It should be borne in mind that several challenges—associated with selection of suitable geometric and radiometric correction equations for large spatial extents, fitting classification methods when the number of classes is large, appropriate interpolation techniques for creating the intensity raster, range ambiguities during data collection, intensity heterogeneity and energy loss (primarily caused due to narrow scan angle), and unidentified influences of laser beam incident angles and illumination of the target material on LiDAR intensity data–already exist, and these issues need to be addressed for fruitful exploitation of multispectral data [35–37, 69]. In addition, the towering cost associated with multispectral LiDAR sensors limits its applicability. Therefore, before acquiring the multispectral LiDAR, the purpose and agenda should be well defined to ensure that the results will justify the investment and meet expectations. In this regard, we recommend the use of multispectral LiDAR on a case-by-case basis, and operations such as land cover classification should be considered as a byproduct or one of the multiple objectives, while using multispectral LiDAR for optimizing the amount spent on data acquisition.

Conclusion In this study, we assessed the capability of a cutting-edge LiDAR sensor, Multispectral Optech Titan, combined with advanced modelling derivatives, for classifying land cover in an urban area. Specifically, we only considered the raster products from LiDAR data (i.e., intensity data and nDSM data), and extracted the MP features from intensity data. A novel hierarchical morphological profiles feature was proposed to extract spatial features of multispectral LiDAR intensity data while maintaining vertical structure information. Results show that the MP feature is useful for providing spatially consistent land cover classification in urban areas. In addition, the proposed HMP feature was found to work best among different features for the multispectral LiDAR data by making use of height information. We obtained an overall accuracy of 93.28% for land cover classification of four classes in the urban area from our best tested model (IMEAN+PseudoNDVI+nDSM+HMP). The results from ourstudy demonstrate that classification results can be greatly improved by extracting spatial features from three-band LiDAR intensity composite images. Future studies could further exploit spectral-spatial classification methods applied to multispectral LiDAR data, and possibly directly classify the point cloud data (i.e., considering geometrical features), which poses new challenges for feature extraction methods.

Acknowledgments This work is jointly supported by the National Natural Science Foundation of China (Grant No. 41401421, 41701397) and the Youth Foundation of President of Institute of Remote Sensing and Digital Earth, Chinese Academy of Sciences (Y6SJ2400CX). We thank Paul LaRocque, Nicola Scotto and Teledyne Optech (http://www.teledyneoptech.com/) team for providing us the sample of the multispectral LiDAR data.