On the Performance of Temporal Stacking and Vegetation Indices for Detection and Estimation of Tobacco Crop

Machine learning in association with remote sensing has assisted agricultural specialists in monitoring, classification and yield estimation of crops. Tobacco is a major taxable crop of Pakistan, however the existing traditional methods for its monitoring and yield estimation are not only expensive and time consuming but also have limitations in terms of accuracy of collected data by a large number of diverse human surveyors. Due to the existence of such loopholes in the employed mechanism for tobacco crop monitoring and yield estimation, its illicit growth and distribution is on the rise. In this paper we have established a sophisticated machine learning mechanism for tobacco crop estimation using temporally stacked sentinel-2 satellite’s data of Pakistan. Instead of the conventional approach of using single remotely sensed imagery for the target crop classification, we propose a machine learning based classification algorithm while keeping in view the phonological cycle of the target tobacco crop. Using the proposed mechanism, the temporal variations within the tobacco crop and its association with the variations of other vegetation is considered to improve the classification performance of the employed machine learning algorithm. Furthermore, the impact of stacking the vegetation indices derived from near infrared and vegetation red edge bands of sentinel-2 with the original sentinel-2 datasets, including Normalized Difference Vegetation Index (NDVI) and Normalized Difference Index 45 (NDI45), on the classification performance of the machine learning mechanism is investigated. Ground Truth data for training of our Artificial Neural Networks classifier, was obtained using indigenously developed survey application “GEOSurvey”. Experiments were conducted using our proposed mechanism while considering various input data setups – including single date imagery, temporally stacked datasets based on phonological cycle of tobacco crop and different combinations of NDVI and NDI45 stacking. Our proposed experimental setup consisting of temporally stacked imagery along with NDVI stacking results in the best classification performance of 95.81% with reference to the single date imagery stacked with NDVI and NDI45, with performance gain of 07.32%.


I. INTRODUCTION
Tobacco is a major revenue generating crop of Khyber Pakhtunkhwa (KP), province of Pakistan. However, the The associate editor coordinating the review of this manuscript and approving it for publication was Yakoub Bazi . government agencies are facing challenges in tobacco crop monitoring and yield estimation, which results in its illicit growth and distribution [1]. Both the national and international cigarette manufacturing companies have very high demand for the tobacco produced within Pakistan. It has been estimated that 21% of the country's Gross Domestic  [8]. The graph represents an increase in trade of illegal cigarettes over the months. The month of December shows trade of 40.6% of unregistered cigarettes. Furthermore a percentage of 85 is due to locally unlicensed tobacco while 15% of it comes due to the result of smuggling [8].
Product (GDP) comes from agriculture [1], [2]. Pakistan produces 100,000 tons of tobacco annually and over 1.2 million people are associated with the tobacco industry [3]. According to an annual report by Pakistan Tobacco Board(PTB), 95% of Flue Cured Virginia (FCV) tobacco produced in the country is harvested in KPK [4]. Despite the fact that this crop is a massive revenue generator for a developing country, the growth and selling of illicit tobacco is on the rise. According to statistics of 2016 by Pakistan Tobacco Company (PTC), 2 out of 5 cigarettes sold in Pakistan are illicit, which in turn incurs a loss of PKR 20 Billion per year [4]. [ Figure 1 [4]] illustrates the market share of illicit tobacco in Pakistan for the year 2016. The graph shows a continuous increase of illicit cigarette trade over the months of 2016. A total of 40.6% illicit tobacco trade is observed for the month of December -2016, in which 85% is due to local DNP(Duty Not Paid) tobacco, while 15% of it is due to smuggling. Despite the fact that tobacco crop is a major contributor to the GDP of Pakistan, authorities in Pakistan are facing challenges in its accurate detection and estimation. Unfortunately, no remote sensing system for tobacco crop monitoring and estimation is established in Pakistan. The existing mechanism employed for tobacco crop estimation and monitoring is based on physical ground surveys and manual measurements using measurement taps or number of footsteps. The existing mechanism of tobacco crop monitoring is indeed very expensive due to the involvement of a large number of human surveyors and error prone.
Currently, remote sensing has been phenomenal in the detection and estimation of different types of vegetation, including their health and status of crop maturity. Many commercial satellites have been launched with a better spectral, spatial and temporal resolution. There exists satellites specifically designed for agricultural monitoring, including drought measurement, vegetation detection and collection of crop health statistics. These satellites include GeoEye-1, Ikonos, QuickBird, Advanced Land Observation System (ALOS), FORMOSAT-2, GeoFen, SPOT 6 and SPOT 7, World View, and many others [5]- [7]. However, for developing countries like Pakistan, the affordability of retrieving such super resolution satellites data is not a cost effective solution to address its agricultural requirements. Fortunately, there are also a number of remote sensing satellites that provide free open remote sensing data for its beneficial utilization. Few of the existing open data satellites are; Landsat (launched by United States Geologi-cal Survey, NASA), MODIS (Moderate Resolution Imaging Spectroradiometer, launched by NASA) and Sentinel (Eu-ropean Space Agency's satellite series). Each satellite has its own specifications and their specifications vary in terms of their spatial resolution, temporal resolution, spectral resolution, revisit time, number of channels, and wavelength of its different bands. Different Channels or bands have different characteristics, designed by engineers for task specific purposes.
In [9] Weiss et al. presented that plant traits and agronomical variables can be estimated by remote sensing and deterministic and empirical approach was adopted for their retrieval. Furthermore in [9] the author provided a detailed description of the recent advancements applied to remote sensing in the field of agriculture. These agronomic traits can be related to various aspects of the crop -e.g. crop type (typological), crop phenology (biological), leaf nitrogen content (chemical), leaf inclination (structural), plant density (geometrical) or crop canopy temperature or soil moisture (physical). Weiss et al. in [9] declared the need of remote sensing technology for agricultural analysis as inevitable and non-destructive. Sonobe et al. in [10] compared indices derived from sentinel-2 for land cover and land-use classification. Furthermore, the performance of Support Vector Machines(SVM) and Random Forest(RF) classifiers have been analysed for land-cover classification. It has been concluded in [10] that the use of vegetational indices can assist in identifying specific crop types. Furthermore, Belgiu and Csillik in [11] Utilized Time Weighted Dynamic Time Warping (TWDTW), for pixel based and object based classification approaches while utilizing RF classifier and their performance is compared. Vuolo et al. [12] examined an agricultural region in Austria and considered 9 crop types within the time span of two years. They assess the impact of multi-temporal information in their dataset. Results showed that the crop type classification accuracy increases with addition of multi-temporal information and similar trends were observed for the year 2016 and 2017.
A very beneficial milestone in the paradigms of open data was achieved in 2015 with the launch of European Space Agency's remote sensing satellite -i.e. Sentinel-2A. Programmed to operate in a couple, its twin satellite Sentinel-2B was launched in 2017. With 10m (Red,Green, Blue, NIR), 20m (Vegetation Red Edge 1, VegetationRed Edge 2, Vegetation Red Edge 3, Vegetation Red Edge4, SWIR-1 and SWIR-2), and 60m(Water Vapor, SWIR-Cirrus and Coastal Aerosol) spatial, 5 days temporal, and 13 bands spectral resolution.
In our proposed research work, we utilize temporal data of Sentinel-2A and 2B for the month of May, 2019 and performed our experiments in collaboration with PTB (Pakistan Tobacco Board). A model based on Artificial Neural Networks (ANN) was developed for tobacco crop detection. Ground Truth data for training of our proposed Artificial Neural Networks classifier, was obtained using indigenously developed survey application ''GEOSurvey''. Various data setups including single date imagery, temporally stacked imagery while considering the phonological cycle of tobacco crop and combination of various vegetation indices, including Normalized Difference Vegetation Index (NDVI) and Normalized Difference Index 45 (NDI45), were considered, in order to investigate its impact on tobacco crop classification performance.
The rest of the paper is organized as follows; Section I provides an overview to the utility of tobacco crop for the selected pilot region of Pakistan along with presentation of state-of-the-art research work in the area of remote sensing and agriculture monitoring and estimation. In section II detailed description is provided about the selected pilot region of Pakistan and penology of tobacco crop and other similar seasonal side crops is presented in order to envision the scope of the target classification work. Furthermore, details about the adopted ANN based classification strategy is presented in (3). Section IV provide detailed information about our adopted strategy for ground truth data collection and pre-processing, followed by description of our adopted experimental setup. Section V presents the findings of our experimentation work. Finally, this paper is concluded with synthesis of the key points of the paper in section VI.

II. RELATED WORK
Tobacco is one of the important commercial crops in several countries including Pakistan, India, China and the United States, with China being the largest producer of this crop [13]. In developing countries the estimation of crops is done manually [14]. Fan et al. used a technique in which tobacco plants are detected and counted using unmanned aerial vehicles (UAV) images [15]. An algorithm has been proposed by the author that is based on the concept of deep neural network. This algorithm consists of three stages mainly tobacco regions extraction, classification and postprocessing. The extraction of tobacco regions takes place through watershed segmentation. The region can have either tobacco or non-tobacco plant. For classification of tobacco plant regions, the training of convolution neural network [16] was conducted. Furthermore, non-tobacco plant regions are removed through post-processing. This is tested on a UAV image data set which includes fourteen tobacco regions. In 2002 Svotwa et al. proposed the use of different vegetation indices for tobacco yield estimation through remote sensing [17]. Keeping in view the importance of this crop new paradigms of research in the detection, classification and field estimation are being carried out today. Zhu et al. perfomed this task using RGB images taken using a UAV [13]. They have conducted experiments on tobacco field using only RGB images, which resulted in the overall accuracy of 95.93%. With the increase in temporal resolution of satellites, temporal analysis and classification is the new technique for identification and estimation of different crops. This can be achieved while knowing the ground truth data, which is crop phonological cycles [18]. Palchowdhuri et al. use random forest classfier [19] different crop types using multi temporal images of world view [20] and Sentinel-2 [18] satellites [21]. The experiment achieved a result of 91% of overall accuracy.

A. AREA OF STUDY
District Swabi is considered to be the main source of tobacco crop. The area we chose for experimentation is Yar Hussain, which is a tehsil of district Swabi. The land of Swabi is mostly plain and agricultural, irrigated through Kabul-Chitral Basin. The soil of Swabi is mostly loamy, and fit for vegetation. The area comprise of a total of 20347.25 Hectares. Figure 2 shows the map for locality of Yar Hussain. There are two types of tobacco plant harvested in Swabi District; White Patta (WP) and Flue Cured Virginia (FCV). Main crops consists of Tobacco, Wheat and a diverse number of vegetation. In Pakistan tobacco and wheat has almost same phonological period. unlike wheat, tobacco seeds are first nurtured in nurseries from the start of November to the end of January. After the nursery period of the plant, cultivation season starts from February and ends in March. Tobacco crop reaches its mature stat in the month of May, while the harvesting season starts in June. Wheat has the same plantation time but no nurseries are used for their nurturing and protection from cold weather conditions. Wheat is also harvested in May and June. Other vegetation such as corn, forage, garlic, ladyfinger, beans, potato, sharsham, shotal, tomato and water melon does not have any specific phonological cycle and are planted randomly. The distribution of fields in the area is random.

A. CLASSIFICATION USING NEURAL NETWORKS
Artificial Neural Networks (ANN), which is inspired from the biological brain, has a motivation to do tasks similar to the human brain. Human brain has various sections to perform complex tasks [22]. Artificial neurons are the building blocks of an ANN. As shown in Figure 3, at least three layers are required to implement an ANN [23]. The first layer is input layer that distributes input to middle layer. The middle layer is hidden that consists of a number of Processing Elements (PEs), used to solve the problem. The number of PEs depends on the complexity of a problem. The last layer is called the output layer, used to generate output on the bases of input parameters.
ANNs have revolutionised machine learning which include pattern recognition [24], classification [25], enhancement [26], analysis [27], estimation and prediction [28]. ANN is   referred as non-parametric because there are no evident sets of parameters that require estimation. Generally, each ANN is trained first, and compared to incoming pixels, by checking on which side of the linear separating surface they lie [28]. A depiction of a typical artificial neuron is shown in Figure 3. The necessary processing element model is based on input from the previous layer and weights. An individual neuron or processing element takes a vector of inputs x = x 1 , x 2 , . . . , x N . Whereas w is connection weight and θ is a bias.

B. FEED FORWARD NEURAL NETWORKS (FFNN)
The feed forward neural network(FFNN) is the simplest type of ANN [29]. FFNN uses a supervised classification scheme. An FFNN is also known as multi-layer perceptron (MLP) [30]. The goal and function of FFNN is to approximate for a classifier. Learning occurs by adjusting the weights in the node to minimize the difference between the output node activation and the input. Approximate mapping of an input pixel into given set of categories (based on the calculated sum of the products of the weights) and the inputs is calculated at each node. Each pixel is processed to make a class selection for each pixel. The network has to learn weights and biases so that the output from the network correctly classifies the input data. Feed-forward networks have the following characteristics: • PEs are arranged in layers, input data travels via first (input) layer and the output layer producing outputs. The intermediate layers have no connection with the external world, and hence are called hidden layers.
• In this network, the information moves in only one direction, namely forward, from the input nodes, through the hidden node(s) and to the output nodes.
• FFNN is acyclic which means that within the network there are no feedback connections among neurons w x to send back the information [31]. There are connections, where a PE h 1 is connected to input x 1 , x 2 and x 3 and h 2 is connected to input x 1 , x 2 and x 3 given in Figure 3. So here all the weights are represented in matrix below (1). The PE computes matrix product of the hidden with those weights plus its own bias then passes to the activation function as well, matrix computation is given as; The mathematics behind the hidden layer is the weighted sum of each input multiplied by the column of input given together as H i (2) (3). The mathematical computation model for neuron is given in Figure 3 which is the computation that will take place on each artificial neuron. All this computation is performed to figure out how to tune the system to match the output to what we think it should be. Small alteration in weight causes significant corresponding change in the output from neural network, this property will make learning possible. But as it occurs a small alteration in weights leads to a big change in the output [15].

C. PARAMETER SETTING FOR FFNN
It is necessary to make several key decisions beforehand. First, the number of layers to use must be chosen. Usually, a three layer network is sufficient, with the purpose of the first layer as an input layer. The inner layers acts as a hidden layer while the output layer is where final decision of a class is made [32]. The next choice relates to the number of elements in each layer. The input layer will generally be given as many nodes as there are components (features) in the pixel vectors. After extensive testing of ANN with a number of changing parameters we found that by increasing the number of layers the classification results in over fitting. Thus we only use a 3 layer network.
The parameters setting of ANN are chosen experimentally such as; training threshold contribution, training rate, momentum rate, training RMS exit criteria, number of hidden layers and training iteration. Algorithm for FFNN is provided in 1.
1) Training Threshold: Training threshold(TT) contribution value lies between 0 and 1 [32]. TT value adjusts the weight of the node. A low value possibly result in inaccurate outputs. We configured TT = 0.7. 2) Training Rate: How small our steps are towards the targeted output? Lowering the training rate is the key to get maximum output.
Step size is any value between 0 to 1. A higher training rate results in computational speedup but with the possibility of unsought results. So, for training accuracy a value of 0.4 is chosen after experimentation [33]. 3) Training Momentum: An aggregate of the small steps taken in the random direction point towards minimum loss of the classifier. This is achieved by keeping the running average of the gradients instead of the direction of current batch of the data. Training momentum can take a value between 0 to 1. It should always be greater than 0 to achieve higher training degree value without oscillating [33]. We selected TM = 0.5 for imagery under study. 4) Training RMS Exit Criterion: RMS error value at which the training should stop. Regardless of whether the number of iteration cycles has not been satisfied [34].

5) Number of Hidden Layers:
The number of hidden layers increases complexity and helps in generating projected results. We chose them 3 [35]. A ground survey to collect ground truth information has been performed in the area of interest. The survey has been performed by visiting the area. Traditional tools used for these types of surveys include global positioning system (GPS) devices, which are expensive and their operation time consuming. In order to overcome these discrepancies we developed our own android based system (Geo Survey) for ground truth collection. With the capabilities to chose the method of survey, the user can have multiple choices to conduct the operation Figure 4(a). Figure 4(b) shows a polygon drawn by choosing the tapping option from the main menu, while Figure 4(c) represents viewing of the survey conducted and polygons drawn. In Figure 4(d) the data viewing capability of the application is shown. The application is freely available on google app store https://play.google.com/store/apps/details?id=com. ncbc.survey.gisGeo Survey.
The developed application is native, which is using JAVA programming language. The survey data is being saved in google firebase real time database. Data from firebase is being downloaded in JSON format, and converted into KML using indigenous python scripting. The database used for the storage of data is MySQL. At the end KML is converted into shapefiles using ARCGIS. With the choice of retrieving a polygon by encircling or by selecting different points interactively, our survey application proved to be cost effective and time efficient, as compared to other traditional methods.

B. PRE-PROCESSING SENTINEL-2 DATA
The retrieval of Sentinel-2 level 2A leads us to pre-processing step, described in Figure 5. Level 2A data is already pre-processed for atmospheric and radiometric correction. For further operations we perform the following steps. During further pre-processing, extraction of area of interest using shapefiles, calculation of vegetation indices (NDI45, NDVI) and layer stacking take place. Where; and

C. EXTRACTION OF REGION OF INTEREST
For the extraction of region of interest, we use shape files retrieved from https://data.world/ocha-fiss/a64d1ff2-7158-48c7-887d-6af69ce21906. As shown in the Figure 2  name the union of samples 'Other Vegetation', to minimize the bias in classification. FCV and WP are combined into one class, and Wheat is kept as an independent class of its own. The reason for keeping wheat as a separate class is because of having this crop in abundance in the region. Training Data in terms of number of polygons of 4 classes are given in Table 2. These polygons were collected using ground surveying mechanism outlined at Section IV, A. These polygons are in the format of ESRI Shapefiles with a smallest polygon representing an area of25m 2 , which are later used to train ANN classifier.

E. EXPERIMENTAL SETUP
The experiments on temporal stacked data and 26 th May 2019 scene have been carried out in a manner described in Table 3.

1) TEMPORAL STACKING
We perform our experiments in a way, such that 4 preprocessed scenes of May, 2019 are temporally stacked with each other, which also include their NDVI and NDI45. The reason for layer stacking of multiple Sentinel-2 scenes can be found by the fact of different sizes of tobacco crops during its growth season based on its reflectance feature in different channels. As the Tobacco crop reaches its maturity with time, its reflectacne increases as depicted in Figure 6(b). This becomes the basis of layer stacking multiple Sentinel-2 scenes, helping identify different levels of Tobacco, Wheat and Other Vegetation. This results in retrieving a large number of training data sample space.

2) STANDALONE SCENE OF 26 th MAY, 2019
Standalone scene of 26 th May, 2019 has been used for classification for the purpose of comparing classifier results with temporally stacked data. We perform classification on a single date scene, which is 26 th May, 2019, because of its high vegetation reflectance, hence showing the peak period of tobacco crop ( Figure 6)(a). Figure 6(a) shows GTD means of the month of May where it can be seen that tobacco reflectance in Near Infrared (NIR), Vegetation Red Edge 1, 2 and 3 bands increases with the advancement in dates. The graph of 26 th May, 2019 shows high reflectance for tobacco.

3) LAYER STACKING OF NDVI AND NDI45
Vegetation indexes which are; NDVI and NDI45 [36] were stacked with the data sets in a manner described in table 3. There are a number of Vegetation indexes used in the studies of remote sensing. Based on the VOLUME 8, 2020   [36]. The use of Normalized Difference Vegetation Index (NDVI) is wide in the field of remote sensing. This vegetation index monitors healthy and non-healthy vegetation. With healthy vegetation, absorbing the red channel while reflecting the Near Infrared wavelength. The difference between the two bands identify healthy and non healthy vegetation areas. With each vegetation type having different NDVI index, the NDVI provides additional support to classifiers in classification. More over the NDVI composes a measurement for the photosynthetic activity and is strongly in correlation with density and vitality of the vegetation. On the other hand NDI45 is more linear with less saturation at higher values than the NDVI.
The region of rapid change in reflectance of vegetation in the near infrared range of the electromagnetic spectrum is referred to as Red Edge (RE) position. Chlorophyll contained in vegetation absorbs most of the light in the visible part of the spectrum but becomes almost transparent at wavelengths greater than 700 nm. Sentinel 2 has band 5, designated at wavelength of 704.1 nm, for the measuring of the Red Edge position of the vegetation Table 1. With rapid change in occurrence of vegetation using the Normalised Difference Index (NDI) a more robust and accurate analysis of healthy and non healthy vegetation can be achieved. The prominent spectral feature of vegetation located between the red absorption maximum and high reflectance in the NIR is RE. Quantification of the RE is often achieved through calculation of the red-edge position (REP) which is recognised as the point of maximum slope along the RE and has been argued to provide enhanced estimates of LCC and canopy chlorophyll content [37]. The reliability of a classifiers cannot be distinguished upon a single parameter like 'overall accuracy'. In order to validate a proposed model its producer accuracy, user accuracy and the overall accuracy must be kept in review [38]. The validation criteria for the classifier has been set according to following parameters; 1) Omission Error: This error points to ground truth data that has been left out in the classification image. This type of error is also called type I error of a classifier.
3) Producers Accuracy: Producers accuracy responds to complement of the error of omission. This is the number of correctly classified samples in the ground truth data, for a specific class. This type of accuracy

V. RESULTS AND DISCUSSION
Experiments using 8 setups explained in section IV-E were conducted using neural networks classifier with 1000 iterations and 3 hidden layers. The training RMS (Root Mean Square) exit criteria has been set to be 0.1 throughout the experimentation. Training threshold and training rate has been set to be 0.0 and 0.2 respectively. The results have been calculated in terms of overall accuracy, producers accuracy and users accuracy. Our findings suggest that the temporal stacking of one month Sentinel-2 level 2A data provides the best results along with NDI45 and NDVI indices (setup 6, setup 7 and setup 8) with above 95% of overall accuracy (Table 6, Table 7, and Table 8). Degraded results having 88% of overall accuracy (Table 4) have been recorded in standalone data of 26 th May, 2019 (setup 1). More over with experiments conducted on a single date imagery with only NDVI stacked (setup 4), the accuracy increases to 90.4592%. Classification map with 4 classes is shown in Figure 7. The producer and user accuracy for Tobacco crop are as follows;

A. PRODUCER AND USER ACCURACY FOR STANDALONE 26 th May 2019 SCENE
The recorded producer accuracy for tobacco by using standalone 26 th May 2019 scene without any VI stacked (setup 1) is 97.15% (Table 4), while the observed user accuracy is 90.75%. Standalone scene of 26 th May 2019 stacked along with NDI45 and NDVI (setup 2) provides producer's accuracy of 96.96% and 92.29% of user accuracy (Table 5).    Similarly producer and user accuracy of 96.14% and 93.23%, respectively has been observed for standalone image staked with only NDI45 (setup 3) ( Table 6). Standalone image stacked with NDVI (setup 4) provides 95.8% and 93.89% of producer and user accuracy respectively (Table 7). In the similar manner Table 9, Table 10, Table 11 and Table 8       Much improved and promising results have been observed by experimenting upon setup (5-8) described in section IV-E. This is due to the temporal staking of 4 retrieved scenes of sentinel-2 data for the month of May. A producer accuracy of 96.25% and user accuracy of 97.38% is recorded in experimenting on temporally stacked data without any VI stacked (setup 5) ( Table 8). For setup 6 (With both NDVI and NDI45 stacked), the producer and user accuracy have been recorded as 98.91% and 94.31% respectively (Table 9). 99.17% of producer and 93.77% of user accuracy respectively are recorded for stacked NDI45 (setup 7), (Table 10). Stacking with NDVI of each scene with May 2019 (setup 8), a producer accuracy of 98.99% and user accuracy of 94.89% is observed (Table 11). Finally To comprehend this slight change in overall, producer and user accuracy, the importance of temporal features should be kept in focus. It has been observed that the best producer accuracy of tobacco is 99.17% in setup 7, while a producer accuracy of 98.99% has been recorded in setup 8, which shows a decrease of 0.8%. Similarly the detected producer accuracy for both the VIs stacked with sentinel-2 data (setup 6) shows a producer accuracy of 98.91%, which also shows a decrease in the producer accuracy by 0.26% as compared to setup 7. During the stacking process vegetation indices were calculated separately for each of the four temporally different scenes of the pilot region and classification decision was carried out based on its cumulative impact. As the response of the VIs is dependent on the availability of chlorophyll contents at the underlying region, therefore its response will be high for all the vegetation types including Tobacco Crop. It has been noticed that employment of such multiple considerations of the VIs for temporally stacked imagery will have impact on improving the accuracy of vegetation classification, in general, but it may result in creating confusion in observation of the underlying crop types. As a result an increase in error of omission will result in increased error of commission of some other vegetation within the classification process. More specifically, while using temporally stacked imagery with both NDVI and NDI45 (setup 7) the producer accuracy of tobacco crop decreases which result in increase in error of omission of the tobacco crop due to high combined response of both considered VIs, which results in match of the spectral signature of tobacco in few cases with other considered vegetation category. Therefore as a result the user accuracy of the other vegetation category decreases as a result of increase in error of commission of other vegetation category -i.e. classification of few cases of Tobacco crop as other vegetation.
Maximum overall accuracy improvement in both single scene scenario and multi-scene stacked scenario is observed while employing NDVI stacking. It can clearly be observed that by employing multi-scene stacking with NDVI results in the best overall accuracy improvement of around 5.98% (for Tobacco only) with reference to single scene without NDVI and NDI45 stacking. Furthermore, the proposed multi-scene scenario with NDVI stacking has overall accuracy of 95.8114% (for all categories), with 7.32% accuracy improvement with reference to the bench marker scheme of single scene imagery without NDI45 and NDVI stacking.

VI. CONCLUSION AND FUTURE WORK
In this study we focused on the application of machine learning and remote sensing to tobacco crop estimation in Pakistan. The traditional methods for monitoring and estimation of tobacco crop are error prone and require a large number of human resources. Since tobacco is one of the most profitable crops for government in terms of tax revenue, there is an urgent interest in exploring more accurate and cost-effective methods for its monitoring and estimation. Application of Machine learning techniques on open data obtained through satellites is one of the best options for monitoring and estimation of tobacco crop, in terms of results accuracy and cost effectiveness. In this study, various techniques such as single date imagery, temporally stacked imagery of Sentinel-2 with reference to the phonological cycle of tobacco, along with each image's NDVI and NDI45 were experimented. The training data for classifier in the form of pixels was acquired through ground survey using our indigenously developed Android application ''GEOSurvey''. Experiments on data sets of temporally stacked images outperformed those on single date data set. More specifically, it is observed that through the employment of neural network classifier on bench-marked single date imagery, an overall accuracy of 88.49% was achieved, which was improved further to 90.45% through the employment of NDVI stacking. In comparison with bench marker single date imagery classification, our experiments on temporally stacked imagery with NDVI using neural network classifier shows improved overall accuracy of 95.81% with 07.32% accuracy improvement with reference to the bench marker scheme. In our future work, our goal is to utilize deep learning based classification mechanism for tobacco crop classification and to compare its performance with our existing model. He has supervised many research projects at undergraduate and graduate levels. His research interests include many aspects of the physical and MAC layers-based cooperative diversity wireless communications protocols, underwater cooperative diversity protocols, and half duplex and full duplex underwater communication protocols. He has authored many journal articles (SCIE and ESCI) and conferences papers. He has received research grants from HEC Pakistan as well as from Higher Education Malaysia. He is a member of Pakistan Engineering Counsel (PEC). He is also a Technical Member of National Curriculum Review Committee (NCRC) Higher Education Commission (HEC), Pakistan, and a Technical Member of Academic Counsel (AC) and Board of Study (BOS) Committee Qurtuba University of Science and IT. He has also been serving as a designated reviewer for many reputed international journals.
ZAHID WADUD received the bachelor's and master's degrees in electrical engineering from the University of Engineering and Technology at Peshawar, Pakistan, in 1999 and 2003, respectively, and the Ph.D. degree from the Capital University of Science and Technology, Islamabad, Pakistan, with the thesis entitled, Energy balancing with sink mobility in the design of underwater routing protocols. He is currently working as an Assistant Professor with the Department of Computer Systems Engineering, University of Engineering and Technology Peshawar. His research interests include wireless sensor networks, energy efficient networks and subsystems, mathematical modeling of wireless channels, embedded systems, and sensors interface. He has published 13 peerreviewed articles in highly ranked journals.
MUHAMMAD ZEESHAN received the B.Sc. degree in computer systems engineering from the University of Engineering and Technology (UET) Peshawar. He is currently pursuing the M.Sc. degree in computer systems engineering from UET Peshawar. He has around six years of experience as developer in web development and mobile app development. Based on his experience, he is also leading development work of web and mobile apps at National Center for Big Data and Cloud Computing (NCBC), UET Peshawar. SUHAIL YOUSAF graduated from Quaid-e-Azam University Islamabad, Pakistan in 2003. He received the M.S. leading to Ph.D. degree from the Computer Systems Group, Vrije Universiteit Amsterdam, The Netherlands. He is currently serving as an Assistant Professor with the Department of Computer Science and IT, University of Engineering Technology at Peshawar, Pakistan. There he worked under the supervision of Prof. Maarten van Steen on the design of a cosmic-ray detector as a wireless distributed system. His research interests include high-performance computing (HPC), development of efficient algorithms for processing of large dynamic graphs and the Internet of things (IoT). Apart from M.S. students, he is currently supervising two Ph.D. scholars in these areas. He is also a CoPI with the Intelligent Systems Design Lab, University of Engineering and Technology at Peshawar, Peshawar, under the National Center of Artificial Intelligence, Pakistan. He was awarded a fully-funded scholarship by the government of Pakistan for his M.S. leading to Ph.D. degree.
ABDUL BASEER QAZI (Member, IEEE) received the M.S. degree in information and communication systems from the University of Technology, Hamburg, Germany, and the Ph.D. degree from the UNU-MERIT, University of Maastricht, The Netherlands. He is currently a Senior Assistant Professor with the Department of Software Engineering, Bahria University, Islamabad, Pakistan. Prior to this, he was an Assistant Professor with the Capital University of Science and Technology, Islamabad. He also worked as an Assistant Professor with CECOS University, Peshawar, Pakistan. His industrial experience includes working for four Fortune 500 companies, IBM, Siemens, Philips Medical Systems, and NXP Semiconductors in Germany and The Netherlands. His research interests include wireless networks, antennas, technology policy, innovation, and entrepreneurship. VOLUME 8, 2020