Seasonal Vegetation Changes in the Malinda Wetland Using Bi-Temporal , Multi-Sensor , Very High Resolution Remote Sensing Data Sets

Small wetlands in East Africa have grown in prominence driven by the unreliable and diminished rains and the increasing population pressure. Due to their size (less than 500 Ha), these wetlands have not been studied extensively using satellite remote sensing approaches. High spatial resolution remote sensing approaches overcome this limitation allowing detailed inventorying and research on such small wetlands. For understanding the seasonal variations in land cover within the Malinda Wetland in Tanzania (350 Ha), two periods were considered, May 2012 coinciding with the wet period (rainy season) and August 2012 coinciding with a fairly rain depressed period (substantially dry but generally cooler season). The wetland was studied using very high spatial resolution orthophotos derived from Unmanned Aerial Vehicle (UAV) photography fused with TerraSAR-X Spotlight mode dual polarized radar data. Using these fused datasets, five main classes were identified that were used to firstly delineate seasonal changes in land use activities and secondly used in determining phenology changes. Combining fuzzy maximum likelihood classification, knowledge classifier and Change Vector Analysis (CVA), land cover classification was undertaken for both seasons. From the results, manifold anthropogenic activities are taking place Corresponding author.


Introduction
Globally, the amount of land cover influenced by anthropogenic activity has grown along with human population and living standards [1].This has put pressure on arable and protected lands, such as wetlands, for agriculture to support these populations.Small wetlands in East Africa have gained increasing prominence as areas for growing crops, driven largely by diminishing and unreliable rainfall on the one hand and growing population with the attendant demand for land and water for food production [2].It is therefore necessary to have a good description of the capacity of these wetlands to sustain these pressures through characterization of the land covers, their variability across seasons and the progression of phenology from one season to another.
Satellite remote sensing provides cost effective means for inventorying and monitoring all types of wetlands but the low wetland accuracies obtained from traditional classification and differentiation methods minimize to a large extent of the utility of such inventories [3]- [5].Wetlands in general present a unique challenge for remote sensing classification due to the complexity and heterogeneity of wetland vegetation making it difficult to arrive at good classification accuracies [6].As part of the Small Wetlands in East Africa (SWEA), 51 small wetlands in Kenya and Tanzania were identified and remote sensing supported approaches used to delineate these wetlands [7].This delineation used optical and microwave remote sensing of medium scale resolution (ENVISAT ASAR and Landsat 5 and 7).Given the coverage of the SWEA project, these products were adequate for providing a synoptic view of the various small wetlands.These identified small wetlands were classified as being of two types; inland valleys and flood plains wetlands [8].In their classification, they found that inland valley swamps dominate the humid mid-hill and highlands whilst floodplains were dominant in the semi-arid and sub-humid lowlands.Inland valley wetlands are generally narrow (0.5 -35 Ha), while floodplain wetlands were larger .
During SWEA-I, aircraft-borne aerial photography was undertaken and this was justifiable since many sites were being imaged [7].In SWEA-II, only Malinda Wetland was mapped and thus aircraft borne aerial photography was impractical given the size of the study area and the costs of a manned photographic campaign.In its place, the Unmanned Aerial Vehicle (UAV) methodology was used to undertake the photography.UAVs are suited for mapping areas of fairly small spatial extent, and coupled with the ability of obtaining very high spatial resolution [9], made this approach the obvious choice for the photographic data collection.The radar data used in this work were from the TerraSAR-X sensor, and the SpotLight mode yielding high spatial resolution imagery.TerraSAR-X radar data has been used for soil moisture determination [10] [11], vegetation biomass quantification [12] among a host of applications.
To obtain maximum benefit from various remote sensing sensors, available image fusion methods have been developed.Image fusion is the combination of two or more different images algorithmically [13].Image fusion aims at the integration of disparate and complimentary information to enhance the information in the image and its interpretation, leading to more accurate information and increased utility.References [13] and [14] provide reviews of various image fusion methods and techniques utilizing multi-sensor, multi-temporal, multi-frequency and multi-resolution data.High spatial resolution imageries such as those obtained from aerial photography or from the very high resolution satellites such as IKONOS provide detailed structural information in addition to the spectral information [15].Given the kind of data utilized in this research, image fusion has been undertaken; combining the UAV derived orthophotos with the TerraSAR-X SpotLight imageries.
In between observation epochs, land cover variations occur and can be classified into two main forms; land cover change where there is a complete change from one cover type to another (e.g.rain-fed and irrigated agriculture), and phenology change occurring within the same land cover class, but caused by the differences in reflectance due to variations arising from the seasonal patterns (e.g.natural vegetation).Land surface phenology (LSP) refers to the seasonal pattern of variation in vegetated land surfaces observed from remote sensing [16].Vegetation phenology (VP) on the other hand refers to the timing of plant processes such as the growing stages (flowering, productivity and leaf senescence) [17].LSP observed from remote sensing is aggregated information at the spatial resolution of the remote sensing sensor.
VP has historically been measured at small dispersed plots that offer high confidence levels but is limited in terms of spatial and temporal coverage [18].Remotely sensed data coupled with statistical modelling have been used to map phenology across spatial extents and fine temporal scales [17].Remote sensing of phenology provides a mechanism to move from plant specific to regional and continental scale studies of phenology.In this work, consideration was made on LSP detected from land change cover analysis.
To help determine the proportion of change corresponding to land cover change and land surface phenology, the Change Vector Analysis (CVA) concept was considered.It has been reported in the literature as a valuable technique for land use/land cover change detection [19].In this approach, a change vector can be described by an angle of change (vector direction) and a magnitude of change from one epoch to another epoch [19].If the gray values for the two images are given as , , , n e e e =  E and ( ) , and n is the number of bands, a change vector is defined as where ΔE includes all the change information between the two epochs for a given pixel and the change magnitude ΔE is computed as The greater ΔE is the higher the possibility of change will be, with the decision on change being on the basis of a threshold.
In the Malinda Wetland and the general Tanga region, the month of May coincides with the period of intense rainfall (wet/rainy season), while August corresponds to the onset of the dry season.Observations in May thus allow the capturing of the rain impacts and the associated activities, while observations in August allow the recording of the response of the wetland and the anthropogenic impacts to the reductions in water availability in the wetland and its periphery.From these discussions, small wetlands exhibit a lot of dynamism in the evolution of land cover and phenology from season to season, and from year to year [20], driven by the fluctuations in the amount and availability of water in the wetland, changes in weather patterns, pressure for agricultural use and their own self-preservation and regeneration efforts.These wetlands have a high potential for agricultural use, among an array of other uses and their characterization i.e. seasonality of land cover, phenology and heterogeneity of land covers, is a key factor to their wise utilization.It is therefore important to adequately delineate, characterize and quantify the changes in land covers, changes in phenology and the attendant areal conversions so that decision makers can be able to see patterns in their utilization.This will in turn inform the decisions on the types of measures to be made to allow the sustained and wise use of these protected resources.
This present research therefore aimed at addressing the following aspects within the Malinda Wetland: (1) Generation of a high resolution Digital Surface Model (DSM).
(2) Develop a methodology that can be used to prepare and process optical and radar datasets for classification and change analysis.
(3) Use bi-temporal, multi sensor data to perform classification of the various land covers present in the wetland.(4) Determine and quantify changes in both land cover and phenology taking place across seasons.

Study Area
Out of the 51 small wetlands identified during SWEA-I, the Malinda Wetland was chosen as the wetland on which more detailed studies would be conducted.It is in the Korogwe District of Tanga Region, Tanzania, within the Pangani flood plain.It lies between latitudes 5˚04'22.8''S-5˚05'52.8''S and longitudes 38˚19'04.8''E-38˚20'38.4''E(354.6 Ha).River Mkomazi, a tributary of Pangani (Ruvuma) River forms the lower boundary of the wetland.To the east of the wetland lie the Usambara West Mountains.
The main economic activity of the surrounding area is rain-fed agriculture with little irrigated agriculture on the wetland's periphery.Figure 1 shows the wetland portion considered in this work, with its areal extents determined by the red boundary.This wetland is one of the floodplain wetlands identified from the wetland classification work done by [7].

Data
Some of the data used was obtained during SWEA Phase I and the rest collected and ordered during the progress of the research.Table 1 shows the various data used, the periods covered and the corresponding spatial resolutions.
We used TerraSAR-X imageries obtained at both horizontal polarization (HH) and vertical polarization (VV).These were obtained using the Spotlight mode (spatial resolution of 1.0 m).The radar data comprised TerraSAR-X datasets which were obtained from the German Aerospace Centre (DLR) under the proposal LAN1393.Radar records target responses due to interaction of microwaves with targets, which response depends on the composition of the said targets [21].Ground truth data was collected by a high resolution hand camera with geo-tagging capabilities.The pieces of information collected in addition to the image were the location and the orientation of the photographer.These pieces of information were useful in the post-classification accuracy evaluations.

Research Workflow
This was organized around three phases: data collection phase; data preparation phase and; data analysis phase.The data collection phase comprises preparation and submission of TerraSAR-X radar remote sensing data proposal, planning and execution of the UAV photography campaign and the concurrent collection of reference (ground truth) data.The data preparation phase had the radar data pre-processing work, orthophotos, mosaic and Digital Surface Model (DSM) generation and the fusion of the optical and radar datasets.The data analysis phase comprised the image classification procedures, change vector analyses and change detection and quantification operations.
Figure 2 presents the preparation and analyses components carried out after data collection.In this part, we broke the work into three parts comprising data preparation part, initial analysis part and advanced analysis part.The final products of the work are a Digital Surface Model and quantification of changes both in cover and phenology.
The UAV used in this research is a parachute carrying a payload of over 8 kg, equipped with a gimbal mounted camera (Figure 3).Based on the objective of generating a high resolution DSM for the study area, a decision was made to assess the best flying height to use to assure a resolution of 0.8 m.To support this effort; two preliminary runs were conducted with two average flying heights of 80 m and 600 m. the flying height of 600 m was found satisfactory and was adopted for both May and August 2012 campaigns.It gave a spatial resolution of 0.2 m which was much better than the desired spatial resolution while the 80 m height gave cm level resolution which would have been overkill given the requirements.During the May 2012 campaign, a total of 2600 photographs were taken, while 6900 photographs were taken during the August 2012 campaign.
Figure 4 show samples of these images.In the sub figures (b) and (d), the areas covered by the higher resolution images (a) and (c) are indicated.The resolution adopted for the orthophotos was 0.8 m.The specifications of the taking camera were: Nikon D300; 18 mm focal length; F 3.48; GPS enabled.

Processing and Analysis
After the UAV field campaigns, the individual photographs were assembled into a bundle-block adjustment software which performed the mosaicking (May 2012-2600 photographs and August 2012-6900 photographs).This was accomplished by using automated internal and external orientation procedures and feature matching   (autocorrelation).Orthophotos and the DSM were generated from this process too.
Unlike optical data, radar data contains a lot of speckle that should be minimized before they can be used [22] [23].A number of speckle filters have been developed in the literature such as the Lee Filter, the Modified Lee Filter, Frost Filter, Kuan Filter and the Gamma Maximum A Posteriori (Gamma MAP) Filter [24] [25].Prior to making a decision on the type of kernel to use in filtering, we considered a variety of combinations as recommended in the literature [19] [26].Using the speckle filtering tool in ERDAS Imagine, two approaches were considered: (1) Use the minimum size of feature of interest to decide on kernel size to use in single operation for example, in the case of a feature about 10m long, a kernel size of 9 × 9 is recommended for TerraSAR-X Spot-Light given its spatial resolution is about 1.0 and (2) use the sequential filtering approach where we start with a 3 × 3 kernel, then apply a 5 × 5 kernel and finally a 7 × 7 kernel.Compared to the Lee Filter and the Modified Lee Filter, the Gamma MAP filter gave better visually appealing results maintaining structure as well as reducing speckle much better in both approaches considered above.This agrees with the findings of [27] and was subsequently used as the speckle filter for the rest of the radar imageries.
Figure 5 shows a small portion of the scene covering a water body and some vegetation in its neighborhood.It is clear from the original image that speckle greatly suppresses the meaningful information.By applying the two approaches above, we get visually comparable results.The sequential filter is slightly better, but the downside of requiring more processing made us settle on the single 9 × 9 kernel for subsequent filtering using the Gamma MAP filter.
In this work, we use Ehlers fusion to combine the horizontally polarized radar with the UAV derived optical orthophotos.The Ehlers fusion is based on an IHS transform coupled with a Fourier domain filtering.The Fourier transform of the optical image and the radar image allows an adaptive filter design in the frequency domain.Using Fast Fourier transform (FFT) techniques, the spatial components to be enhanced or suppressed can be directly accessed.The optical spectrum is filtered with a low pass filter (LP) whereas the radar spectrum is filtered with an inverse high pass filter (HP).After filtering, the images are transformed back into the spatial domain with an inverse FFT and added together to form a fused intensity component from the optical image and the radar image.This new intensity component and the original hue and saturation components of the optical image form a new HIS image.As the last step, an inverse IHS transformation produces a fused RGB image.A complete description of the algorithm can be obtained from [28].
When collecting imagery by aerial means, there are variations in illumination that occur and these present a challenge during interpretation and analysis.To help reduce the effects of varying illumination in the resultant mosaic, spectral ratios were considered to be the surrogate images.Prior to performing image ratioing and image fusion, we performed a histogram match using the May 2012 image as the reference.While their use has been limited to vegetation indices, spectral ratios have been reported as being capable of reducing the effects of variations in illumination conditions [21].

Classification
The hard spectrally based classifiers have been criticized as being unsuited for classification using high resolution datasets [5] [6] [15] largely due to their inability to capture the structures in imagery datasets.One of the strategies of improving the classification skill is to introduce the notion of fuzziness in the assignment of classes from the default crisp assignment [29].Given the heterogeneous nature of the wetland vegetation, the fuzzy maximum likelihood classifier was used.For each pixel under consideration, the assignment was to three classes with a distance parameter giving the level of confidence of the assignment to each class.Fuzzy convolution was then used to perform the final allocation of a pixel to the most suitable class.This progresses as follows: using a multilayer classification and respective distances in feature spaces for each classification layer, the convolution creates a new single output layer by computing a total weighted distances for all classes in the window [30].Classes with lower distance values remain unchanged while larger distance values can change to neighboring values if there is a sufficient number of neighboring pixels with lower distance values.
The initial wetland classes considered were 13 and these were combined to capture the main land uses of interest to SWEA.These class combinations are listed in Table 2.
These five classes were used as input into a knowledge classifier.Additional information fed into the knowledge engineer was the wetland contour.This contour was derived from the DSM (elevation of 346m).Based on this, any predominant wetland categories were reassigned to either human influence or permanent natural vegetation.For covers falling in the wetland contour, the classification results from the fuzzy maximum likelihood classifier were used.The generated knowledge engineer was applied to the classification results.Figure 6 shows the Knowledge engineer hypotheses and rules employed in refining the classification results.

Results and Discussion
Figure 7 shows the image mosaics developed from the UAV photo assemblage.May 2012 was a wetter month and within the wetland near R. Mkomazi it is virtually wetter with water collecting in puddles and ponds.August 2012 was a drier month (not raining) and much of the water has drained, with crops and other vegetation being lush leaved hence the greenness observable across the study area.The DSM captures the variations in the study area and there is a clear demarcation of the actual wetland.This is a surface model capturing the top of canopies and only captures ground conditions for bare or sparsely vegetated regions.There are some raised portions in the wetland, serving as homesteads for some of the inhabitants.There was some challenge in the August image in the form of a shadowed band in the middle of wetland, adding challenge to the utilization of the images.
Figure 8 shows the fused images, which now contains contributions from both the UAV images and the radar images (HH polarization).We have used the HH polarization and combined it with the 3 bands from the UAV.From the fused images some of the variations in cover and phenology are visually enhanced.The impacts of human activities and the reduction in water collection are clearly discernible.
Figure 9 shows the spectral ratios computed for the two periods under study.The indices were calculated according to Equation (3) and were computed using the resultant bands from the fused images.
The color assignment in the figure was according to the indices and was as follows: I 1 (red), I 2 (green) and I 3 (blue).Due to the large range in the intensity values in the 3 bands of the May image due to the presence of water during this period, there is an attendant crispness in the delineation of the various covers present.Considering the August image, there is a smoother distribution in the ratios.This can be attributed to the narrowing of the ranges within bands due to the abundance of various crops and vegetation that are at comparable stages in the growth and development process, corresponding to the high foliage duration.It is clear from the August image, the significance of using ratios to suppress variations in illumination conditions, since this image is more interpretable in the regions shadowed.
Figure 10 shows the results of knowledge classification.Due to the presence of inundated conditions in the wetland from the May image, there is markedly reduced human influence, and open water masses can be identified.In this period the wetland is largely inaccessible and where accessible, it is used for rain-fed rice growing.On the other hand, the August image shows increase in human activities largely due to the draining of the wetland.This arises from the pressure from the local people to utilize the wetland for farming activities.
To assess the quality of the classification exercise, classification accuracy assessment was undertaken.Tables 3-5 show these accuracy statistics.For the assessment, 200 samples were used for both classification images.
The 200 samples collected as reference were used to compare the classification results for each of the ground truth points against the actual land cover.The accuracy assessment then comprises a tally of correctly classified samples against the misclassified ones.Table 3 shows fairly good accuracy in determination of land covers (88.5% overall accuracy and 82.9 % K-hat statistics).Table 4 also gives indication of good accuracies (89.9% overall and 84.3% K-hat statistics).
The CVA carried out in this research utilized the indices presented in Equation (1), specifically I 1 (VI) relating the red and the green bands, and I 3 (WI) relating the blue and the green bands.In the absence of a Near Infrared band, I 1 can be considered to carry similar information to the Normalized Difference Vegetation Index (NDVI) while I 3 carries information on greenness and wetness, which can be used as indicators of the soil surface conditions.Reference [31] used NDVI and Bare Soil Index (BI) as the components in a CVA change detection work.
Figure 11 shows the results of the magnitude and the change direction.A threshold of 0.130 was selected as the magnitude above which this would be deemed to be change.This threshold was arrived at based on ground verified changes observed during the two flight campaigns.Some of the ground truth points had some variation in both cover and phenology and this information were used to validate the choice of the value of 0.130 as the threshold.
From the magnitude image, it can be noted that over most parts of the study area, there were variations in intensities between the two periods, driven by the variations in response due to changes in target conditions and characteristics.The change direction image quantifies the indicators for variations in target responses.Visual inspection indicates that most parts within the wetland experienced drying up while portions in the periphery experienced increased foliage, hence the increased greenness.
Figure 12 shows the results of land cover change and change in phenology analysis.The upper portion of the figure gives the overall changes with insets showing a magnification of the upper portion of the wetland.In this part it can be noted that there has been a conversion to the left part of a portion to permanent natural vegetation but overall there has been a lot of human influence between May and August.This can be attributed to the transition from wet and inundated to drier and accessible with the abating of the rainfall season in the region.This character is also identifiable from the phenology change map.In this case, most of the human influence portions show increase in vegetation due to the growing process of planted crops.In May most parcels have just been freshly planted while in August the crops have grown and substantially and increased vegetation (greener and wetter).
While in the wetland, in May most of the areas are inundated or have high water content while in August, most of the wetland has dried up and thus there is reduced water in the wetland.This makes the wetland accessible.Some of the land parcels inside the wetland that were used for rice farming are now drier and the phonological change is clear.There has been conversion to permanent natural vegetation arising from the fact that during the      rainy spell in May, much of the wetland is inundated and this change in the rain depressed period in August.Thus the background is no longer wet but dry.From analysis of the change detection results, it can be seen that the upper middle portion of the wetland has a lot of human influence.This is largely attributable to the fact that in this area, we have some springs which assure availability of water even during the dry periods and this would prove handy for irrigated farming.The changes in land cover and in phenology are summarized in Table 5.From these statistics, it is clear that much of the land cover changes occurring are conversions due to anthropogenic activities (Human Influence) while much of the conversion to permanent vegetation was due to conversion from inundated state and thus the high conversion rates.Much of the anthropogenic activities in the area include clearing of wetland portions and farming.This was driven by the drying up of the draining of the wetland hence increased accessibility.The lowest conversion is to water masses due to the fact that the trend was from wetter to drier conditions meaning reduced chance of water collecting to form masses.It is also evident that there is a lot of cover change in the upper middle portion of the wetland.This is attributable to the fact that in this region, there are a number of springs which provide a dependable source of water for irrigated farming.It can be seen again that within the wetland, there are phenological changes due to fluctuations in water supply with inundated rice (paddy) fields giving way to irrigated rice growing during the drier period in August.
However, much of the phonological changes are within the periphery due to the changes occasioned by the stage in the growing process of vegetation (crops and other vegetation types).By considering the phenology as a proportion of the respective land cover change a trend seems to emerge, other than the permanent vegetation (both trees & shrubs and the Permanent natural vegetation) the other classes show some correlation with the phenology being about 42% of the corresponding land cover change.
This also indicates the variability of the situation in the wetland with the more stable covers (the permanent) confirming the stable state of the phonological response (lower variability) and the more seasonal covers indicating a corresponding higher variability with respect to the change in phenology within the cover.
Figure 13 gives the schematic on the fluctuations in anthropogenic pressures within and around the wetland with respect to the season.In May, much of the wetland is inundated while the regions in the periphery are not and thus a lot of farming activities are restricted to the periphery of the wetland.However, during the dry season from August, the periphery dries up and most of the crops are maturing, while the wetland becomes accessible for farming activities.There is thus a cyclic behavior in the use of the wetland and its periphery between the dry and the wet seasons.

Conclusions
This paper sets out to perform classification of land use and land covers in the Malinda Wetland, one of the Small Wetlands in Eastern Africa, and determine and quantify changes in land cover combined with changes in phenology.The objectives have been met, and in the process, it has been demonstrated that UAV photography can be used to photograph and produce high resolution orthophotos and DSM in support of mapping needs spanning relatively small areas.By combining radar and optical datasets, land cover classification has been accomplished.We have demonstrated that by considering the extent of the wetland portion, through use of the fuzzy maximum likelihood classifier coupled with the knowledge classifier (akin to the decision tree classifier), there has been improvement in the classification results.Accuracy assessment results for the two epochs are respectively 89.5%; 82.94% and 89.85%; 84.31% for both overall accuracy and K-hat statistics for the May 2012 campaign and August 2012 campaign.
From the application of CVA, we have been able to determine and quantify both land use/land cover changes, and changes in phenology.Most of the land cover conversion taking place has been to the human influence class, reflecting anthropogenic pressures (63.01 Ha).It is also clear that most of the changes in phenology are also occurring within this class (26.6 Ha).This is attributable to the fact that most of the cropping has grown and largely covered the ground (hence the greener and wetter perception) in the phenology analysis.We have also shown the seasonality of the anthropogenic activities between the wetland and its periphery.On the whole, we have presented in addition a methodology that can be used for processing high spatial resolution optical and radar datasets for classification, land cover change and phenology detection and quantification.

Figure 4 .
Figure 4. Quality of individual photographs obtained from UAV. Sets (a) and (c) were captured at a general flying height of 80 m, while (b) and (d) were obtained from a flying height of 600 m.

Figure 5 .
Figure 5. Speckle filtering strategy comparison: (a) Represents the original image portion; (b) The results after applying a sequential filtering approach; (c) After the single large filter corresponding to large sizes identifiable.

Figure 6 .
Figure 6.Knowledge engineer hypotheses and rules.

Figure 7 .
Figure 7. Image mosaics and DSM from the aerial and UAV flight campaigns.

Figure 8 .
Figure 8. Fused images across the two time periods.

Figure 9 .
Figure 9. Spectral ratios computed using the indices from the fused images.

Figure 11 .
Figure 11.Change Vector Analysis magnitude and change direction results.

Figure 12 .
Figure 12.Land cover change and change in phenology.

Figure 13 .
Figure 13.Direction of human activities observed in the study area.

Table 1 .
Datasets used and the periods collected.

Table 3 .
Accuracy assessment for classification results of may 2012.

Table 4 .
Accuracy assessment for classification results of August 2012.

Table 5 .
Land cover change and land surface phenology statistics.