Mapping and Monitoring Small-Scale Mining Activities in Ghana using Sentinel-1 Time Series (2015–2019)

Illegal small-scale mining (galamsey) in South-Western Ghana has grown tremendously in the last decade and caused significant environmental degradation. Excessive cloud cover in the area has limited the use of optical remote sensing data to map and monitor the extent of these activities. This study investigated the use of annual time-series Sentinel-1 data to map and monitor illegal mining activities along major rivers in South-Western Ghana between 2015 and 2019. A change detection approach, based on three time-series features—minimum, mean, maximum—was used to compute a backscatter threshold value suitable to identify/detect mining-induced land cover changes in the study area. Compared to the mean and maximum, the minimum time-series feature (in both VH and VV polarization) was found to be more sensitive to changes in backscattering within the period of investigation. Our approach permitted the detection of new illegal mining areas on an annual basis. A backscatter threshold value of +1.65 dB was found suitable for detecting illegal mining activities in the study area. Application of this threshold revealed illegal mining area extents of 102 km2, 60 km2 and 33 km2 for periods 2015/2016–2016/2017, 2016/2017–2017/2018 and 2017/2018–2018/2019, respectively. The observed decreasing trend in new illegal mining areas suggests that efforts at stopping illegal mining yielded positive results in the period investigated. Despite the advantages of Synthetic Aperture Radar data in monitoring phenomena in cloud-prone areas, our analysis revealed that about 25% of the Sentinel-1 data, mostly acquired in March and October (beginning and end of rainy season respectively), were unusable due to atmospheric effects from high intensity rainfall events. Further investigation in other geographies and climatic regions is needed to ascertain the susceptibility of Sentinel-1 data to atmospheric conditions.


Introduction
Mining (especially surface) is one of the major causes of land and environmental degradation globally [1][2][3]. Environmental impacts such as deforestation, landscape degradation, alteration of stream and river morphology, widespread environmental pollution, siltation of water bodies, biodiversity loss, etc., have been noted to be associated with mining [1,4]. Despite these environmental impacts, mining significantly contributes to the economies of most developing countries. In Guinea and South Africa, for example, Aryee [5] reported that mining contributes 25% and 5.9%, respectively of their Gross Domestic Product (GDP).
The mining sector in Ghana has always been an important contributor to the economic development of the country. Mining is among the top foreign direct investment earners for Ghana [6,7]. From 2011 to 2013, the sector contributed between 8.4 and 9.8% to the country's GDP [8]. Together with cocoa,

Study Area
The study area is located in South-Western Ghana and encompasses an area of 9071 km 2 ( Figure 1). It overlaps three administrative regions (Ashanti, Central and Western South) which are the main regions in which industrially operated gold mines in the country are located [12]. Some of the worst affected municipalities and district assemblies such as Amansie West and Central, Wassa Amenfi East and Central, Wassa West, Upper Denkyira, Obuasi Municipal and Adansi North are located in, or overlap, the study area. Major rivers that drain the area, along which illegal mining takes place include Ankobra, Offin, Oda and Pra.
The largest settlement in the study area is the mining town of Obuasi, which is home to Ghana's largest gold mine AngloGold Ashanti. Two other industrial mines that undertake surface mining are located northeast of Bogoso and west of Ayanfuri ( Figure 1). The study area falls in the moist semi-deciduous forest agro-ecological zone of the country. Average annual temperature is about 26 • C. The agro-ecological zone is characterized by two distinct rainfall seasons: a dry season that spans from December to February and a wet season from March to November. Unlike temperature, precipitation is relatively variable within the study area. For example, from north to south, annual precipitation of prominent towns are 1450 mm (Obuasi), 1531 mm (Dunkwa) and 1682 mm (Bogoso) [28]. The study area is covered by both evergreen and deciduous forests. The latter lose their leaves for a short period of the year, whereupon new ones immediately begin to grow [29]. The largest settlement in the study area is the mining town of Obuasi, which is home to Ghana's largest gold mine AngloGold Ashanti. Two other industrial mines that undertake surface mining are located northeast of Bogoso and west of Ayanfuri ( Figure 1). The study area falls in the moist semideciduous forest agro-ecological zone of the country. Average annual temperature is about 26 °C. The agro-ecological zone is characterized by two distinct rainfall seasons: a dry season that spans from December to February and a wet season from March to November. Unlike temperature, precipitation is relatively variable within the study area. For example, from north to south, annual precipitation of prominent towns are 1450 mm (Obuasi), 1531 mm (Dunkwa) and 1682 mm (Bogoso) [28]. The study area is covered by both evergreen and deciduous forests. The latter lose their leaves for a short period of the year, whereupon new ones immediately begin to grow [29].

Sentinel-1 Time Series
A total of 155 S-1 scenes covering the period 23 July 2015 to 9 April 2019 were downloaded from the Sentinel data hub (scihub.copernicus.eu) as Single-Look-Complex (SLC) products. The scenes were acquired in the Interferometric Wide Swath (IW) mode with VV and VH polarizations and were available in a 12-day interval. All datasets were acquired at 6:25 pm local time in the evening. Although all the datasets lie in the same relative orbit (number 45), there are slight differences in terms of location/extent. For example, for the period 20 November 2015 to 13 May 2017, a more northerly and southerly scene of one day each had to be used to ensure that the study area was completely covered. Consequently, the 155 scenes cover 117 unique days. Scenes acquired at the same data were mosaicked after the processing. However, not all scenes could finally be used in the analysis due to atmospheric effects that reduced the image quality (see Section 3.1). The data were processed using the Sentinel Applications Platform (SNAP) software (version 7.01.) (step.esa.int/main/toolboxes/snap/), Remote Sens. 2020, 12, 911 5 of 26 IDL (version 8.7), and ENVI (version 5.5) (www.harrisgeospatial.com). The processing of the S-1 time series is illustrated in Figure 2 along with the general workflow of the approach. The processing included standard processing steps (i.e., Apply Orbit File, Split, Deburst, Calibration, Multilooking, Terrain Flattening, Boxcar Speckle Filter using a moving window of 5 by 5 pixels, and Terrain Correction). The data were projected to a common gird (World Geodetic System (WGS) 1984, Universal Transverse Mercator (UTM) Zone 30 North) with a pixel spacing of 20 by 20 m, which is not the best possible resolution achievable (approximately 12-15 m), but a compromise considering the speckle effect and the radiometric quality. The terrain corrected gamma nought intensities of the VV and VH channel were used in the subsequent analysis. From the full time series, temporal subsets were generated. These sets were defined as follows: terms of location/extent. For example, for the period 20 November 2015 to 13 May 2017, a more northerly and southerly scene of one day each had to be used to ensure that the study area was completely covered. Consequently, the 155 scenes cover 117 unique days. Scenes acquired at the same data were mosaicked after the processing. However, not all scenes could finally be used in the analysis due to atmospheric effects that reduced the image quality (see Section 3.1). The data were processed using the Sentinel Applications Platform (SNAP) software (version 7.01.) (step.esa.int/main/toolboxes/snap/), IDL (version 8.7), and ENVI (version 5.5) (www.harrisgeospatial.com). The processing of the S-1 time series is illustrated in Figure 2 along with the general workflow of the approach. The processing included standard processing steps (i.e., Apply Orbit File, Split, Deburst, Calibration, Multilooking, Terrain Flattening, Boxcar Speckle Filter using a moving window of 5 by 5 pixels, and Terrain Correction). The data were projected to a common gird (World Geodetic System (WGS) 1984, Universal Transverse Mercator (UTM) Zone 30 North) with a pixel spacing of 20 by 20 meters, which is not the best possible resolution achievable (approximately 12-15 meters), but a compromise considering the speckle effect and the radiometric quality. The terrain corrected gamma nought intensities of the VV and VH channel were used in the subsequent analysis. From the full time series, temporal subsets were generated. These sets were defined as follows:

SRTM Digital Elevation Model
The 30 m resolution digital elevation model (DEM) from the Shuttle Radar Topographic Mission (SRTM) was downloaded from the Earth Explorer portal of the United States Geological Surveys (USGS) (earthexplorer.usgs.gov). The SRTM data was used to conduct hydrological/watershed analysis that can be used to locate activities of illegal miners in the study area in terms of hydrographical position (Figures 1 and 2). The DEM underwent hydrographic modeling in ArcMap 10.6.1. The "Fill" function was applied during hydrographic pre-processing to remove sinks. Flow direction and flow accumulation were subsequently calculated. The channel network was deduced using a watershed size of >4 km 2 . Using this channel network, the stream order of Shreve [30] was calculated for each line segment. This allows us to characterize detected activities in terms of stream order, e.g., if activities are found on small or large tributaries. Finally, the watersheds of the Ankobra and the Pra rivers were approximated by using the "Watershed" function and pour-points placed at the lower southern border of the study area.

Reference Data
Reference data were extracted from two sources including time-series land use/land cover data (LULC), and high-resolution Google Earth Imagery of a smaller area that was strongly affected by mining activity between 2015 and 2017. The LULC data, obtained for the nominal years 2006 and 2015, were derived from regional scale analysis by Gessner et al. [31] and Forkuor et al. (in preparation) based on multi-sensor satellite data. Both sources of reference data were used to identify changed and unchanged targets for the general scattering assessment and subsequent determination of thresholds for detecting mining activities (see Figure 2).
The high-resolution Google Earth images used as source of reference were acquired in February/April 2015 and May/December 2017 ( Figure 3). These images clearly indicated mining activities along a small tributary of the Ankobra River in the western part of the study area, close to the town of Diaso. The reference data covered an area of approx. 19 by 14 km (266 km 2 ) and 104 locations, affected by mining activity, were mapped based on visual image interpretation. On average, these polygons had a size of approximately 3.5 ha, while the smallest changed area that was mapped had a size of 0.5 ha (i.e., smaller areas than this were not considered). It was not possible to compile further meaningful and reliable high-resolution images of larger areas within test site. This was due to cloud cover and/or unsuited acquisition dates, e.g., most recent high-resolution images available in Google Earth date to 2015 and are, therefore, not useful for referencing of changes, as S-1 acquisition started in July 2015.

Change Detection
Change detection technique was used to map the extent of illegal mining along water bodies in the study area based on Sentinel-1 time-series data. The methodology is premised on the fact that analysis of SAR images acquired with the same acquisition geometry (i.e., same relative orbit) will reveal changed areas as having different scattering characteristics over time, while unchanged areas will have very similar scattering characteristics over time. Comparing the intensities between two sets, the variations in intensity of unchanged targets should be low and approximately in the order of +/-1 dB, as reported by References [32,33]. Further, as the establishment of new mining areas usually infers a land cover change from vegetation to bare areas or water, the magnitude of the change in intensity is expected to be fairly high, exceeding the above mentioned 1 dB criteria. It is, therefore, reasonable to attribute significant differences in SAR intensities between two image periods to changes in land cover. In this regard, our methodology assumed that detected land cover changes along major rivers in the study area are induced by illegal mining activities, and the magnitude/extent of the change represents the area under illegal mining. The change detection methodology entailed four steps.

Sensitivity of Sentinel-1 Time-Series to Land Cover Changes
Sentinel-1 image time-series were analyzed to ascertain its suitability to detect changes in land cover in the study area. For each of the four temporal subsets defined in Section 2.2.1 and polarization channels (VV and VH), three time-series features were extracted along the temporal dimension of the

Change Detection
Change detection technique was used to map the extent of illegal mining along water bodies in the study area based on Sentinel-1 time-series data. The methodology is premised on the fact that analysis of SAR images acquired with the same acquisition geometry (i.e., same relative orbit) will reveal changed areas as having different scattering characteristics over time, while unchanged areas will have very similar scattering characteristics over time. Comparing the intensities between two sets, the variations in intensity of unchanged targets should be low and approximately in the order of ±1 dB, as reported by References [32,33]. Further, as the establishment of new mining areas usually infers a land cover change from vegetation to bare areas or water, the magnitude of the change in intensity is expected to be fairly high, exceeding the above mentioned 1 dB criteria. It is, therefore, reasonable to attribute significant differences in SAR intensities between two image periods to changes in land cover. In this regard, our methodology assumed that detected land cover changes along major rivers in the study area are induced by illegal mining activities, and the magnitude/extent of the change represents the area under illegal mining. The change detection methodology entailed four steps.

Sensitivity of Sentinel-1 Time-Series to Land Cover Changes
Sentinel-1 image time-series were analyzed to ascertain its suitability to detect changes in land cover in the study area. For each of the four temporal subsets defined in Section 2.2.1 and polarization Remote Sens. 2020, 12, 911 8 of 26 channels (VV and VH), three time-series features were extracted along the temporal dimension of the time-series. These features were the minimum, mean, and maximum. To evaluate the sensitivity of Sentinel-1 data to land cover changes, the VV and VH scattering characteristics (histograms) of the time series features minimum, mean, and maximum were derived for two broad land cover classes: "water" and "agro-forest". These two classes were selected because illegal mining activities in the study area usually cause a removal of vegetation and the creation of dug-outs or pools of water. Water is thus the dominant land cover after vegetation removal, between which there could be bare soil or vegetation. It is assumed that such changes in land cover (i.e., vegetation to water) will be highlighted in the time series features of S-1, as vegetation and water have different scattering characteristics and the received intensity from vegetation is usually higher than the intensity from water surfaces.
Reference information for the two LULC classes were extracted from existing LULC datasets of Gessner et al. (2015) and Forkuor et al. (in preparation) (see Section 2.3). The overall separability between the two classes offered by the respective time series features (i.e., minimum, mean, maximum) was investigated based on the extracted statistics of the classes and two metrics "distance" (DIST) and "similarity" (SIM) as defined by [34]. Both metrics quantify the distance/similarity between two probability functions (i.e., the histograms of the classes "water" and "agro-forest"). Higher values of DIST and lower values of SIM point to a better separation offered by the feature under investigation, e.g., how well the histograms of both classes separate in the VV mean intensity of set "2015/2016". Therefore, the two metrics DIST and SIM enable a determination of the time series feature(s) that is/are most sensitivity to mining-induced land cover changes. Further, the class statistics indicate the expected magnitude and direction of changes e.g., when a land cover changes from "agro-forest" to "water" occurs.

Determination of Threshold for Detecting Changed and Unchanged Areas
The time series feature(s) (i.e., minimum, maximum, mean) identified as the most sensitive to land cover changes were further analyzed in order to find a threshold that is suited to identify/classify "changed" and "unchanged" areas in the S-1 time series features. First, changed and unchanged areas were mapped from high-resolution Google Earth Imagery. Visual interpretation was used to identify areas that had experience change due to mining activities between 2015 and 2017, as well as unchanged areas. Second, difference images between the time-series features of the temporal subsets "2015/2016" and "2016/2017" (i.e., the difference images of the VV/VH minimum/mean/maximum intensities) were derived. Based on an overlay of the difference images and the mapped "changed" and "unchanged" areas, histograms of the "changed" and "unchanged" areas were derived and analyzed to determine a suitable threshold for detecting mining-induced land cover change between the two temporal subsets (i.e., "2015/2016" and "2016/2017").
Two approaches were tested to determine a threshold. The first is a supervised technique where the threshold is determined based on the intersection of the two histograms (i.e., the intensity value where both "changed" and "unchanged" classes show the same frequency). The second is an unsupervised methodology proposed by Otsu (1979), which iteratively seeks to minimize the interclass variance. The Otsu algorithm was directly applied using the difference images between the time-series features (i.e., no reference information on "changed" and "unchanged" areas is needed to compute the threshold).
Both approaches realize a simple form of binary image classification using a one-dimensional feature space and perform well if bimodal distributions are investigated [35], which was the case here. However, they do not necessarily identify the same threshold, which can provide an insight into the influence of changing thresholds in detecting mining-induced-changes. Additionally, the two metrics DIST and SIM were calculated to ensure that the feature offering the best separation is used in the subsequent step of classifying/mapping changes.

Detection of Illegal Mining Areas
Thresholds determined from the two approaches (i.e., histogram intersection and Otsu) (Section 2.3.2) were applied to difference images computed from the temporal subsets to identify changed and unchanged areas. Difference images were computed for three periods as follows: Period 1: "2015/2016"-"2016/2017", Period 2: "2016/2017"-"2017/2018", and Period 3: "2017/2018"-"2018/2019". The application of the threshold values to the difference maps resulted in multiple change/unchanged binary maps which were subsequently validated using reference data. Prior to this validation, a majority filter with a moving window size of 5 by 5 pixels was applied in the post-classification in order to remove isolated single pixels.

Accuracy Assessment and Analysis
For Period 1 ("2015/2016"-"2016/2017"), the change/no-change maps generated (based on threshold, polarization, time-series feature) were validated using reference data in order to assess the accuracy of the change detection and to indicate the most suited time series feature and polarization. Reference data of changed areas were derived from Google Earth Imagery (see Section 2.3) and formed the basis of this assessment. Producer's and user's accuracies of the changed class, as defined by Congalton and Green [36], were computed. In addition, errors of omission and commission were computed and presented. As several change/no-change maps were produced (i.e., for each polarization, time series feature and threshold), the accuracy assessments provided the basis for selecting the final map for estimating the temporal evolution of illegal mining area extent in the study area. Finally, the spatial pattern and temporal evolution of the changes were analyzed, also in relation to the stream order, which was derived from the SRTM DEM and linked with the classification results via a spatial join in ArcMap (10.6.1).

Atmospheric Effects on Sentinel-1 Data
Our analysis revealed that a significant portion of the investigated S-1 data showed locally-restricted patches of strong local decrease in backscatter intensity in the order >−10 dB compared to the long-term mean intensity. Such areas (approx. 20-50 km 2 in size) of decreased intensity were found over the entire imaging area and did not show a specific spatial distribution or pattern ( Figure 4). Further, these patches are found in single scenes only (i.e., not in consecutive scenes). The shapes of the patches are diffuse and borders are not clearly distinguishable ( Figure 4). A total of 29, out of 117 (ca. 25%), scenes showed this effect and accordingly these scenes were not used in subsequent analyses. While no such patches were found in the 29 scenes acquired in August, in almost two thirds (7 out of 11) of the scenes acquired in March such patches were found. In October, half of the downloaded S-1 images (5 out of 10) were removed from the subsequent analyses. May and

Feature Selection and Threshold Definition
The separability offered by the time series features minimum, mean, and maximum was estimated for both polarization channels for the set "2016/2017" using the "water" and "agro-forest" classes from the two LULC classifications (see Figure 5). The minimum and mean features showed very clear separation of the values of both classes (i.e., water and agro-forest), while maximum

Feature Selection and Threshold Definition
The separability offered by the time series features minimum, mean, and maximum was estimated for both polarization channels for the set "2016/2017" using the "water" and "agro-forest" classes from the two LULC classifications (see Figure 5). The minimum and mean features showed very clear separation of the values of both classes (i.e., water and agro-forest), while maximum features (at VV and VH) showed less separability but rather increased overlap of classes' frequency distributions.
Remote Sens. 2019, 11, x FOR PEER REVIEW 11 of 26 features (at VV and VH) showed less separability but rather increased overlap of classes' frequency distributions. These findings are also supported by the DIST and SIM metrics (Table 1). Lowest SIM and highest DIST values were obtained for the minimum features, followed by the mean features. The maximum features showed higher SIM and lower DIST values. Both polarizations performed similar and values of the metrics were nearly the same for both VV and VH. Therefore, the minimum and mean features of both polarizations offered a better separation of the classes "water" and "agroforest" compared to the maximum features, while there is no clear difference in performance between VV and VH. Consequently, the minimum and mean features were found to be more sensitive to mining-induced changes than the maximum. Table 1. Results of the distance (DIST) and similarity (SIM) metrics of classes "water" and "agroforest" for the investigated time series features minimum (MIN), maximum (MAX), mean (MEAN) and for both polarization channels (VV/VH) (see Figure 6). The two metrics quantify the separation between the frequency distributions of both classes. Higher values in DIST and lower values in SIM point to a better separation offered by the feature under observation. Difference images of the VV/VH minimum/mean intensities between periods "2015/2016" and "2016/2017" were used to derive thresholds for detecting changes in illegal mining areas. Results of this analysis are shown in Figures 6 and Table 2. The high-resolution Google Earth images (see Figure  3), available for a smaller portion of the study area, were analyzed to identify targets ("changed"/"unchanged"), which served as the basis for evaluating the time series features minimum These findings are also supported by the DIST and SIM metrics (Table 1). Lowest SIM and highest DIST values were obtained for the minimum features, followed by the mean features. The maximum features showed higher SIM and lower DIST values. Both polarizations performed similar and values of the metrics were nearly the same for both VV and VH. Therefore, the minimum and mean features of both polarizations offered a better separation of the classes "water" and "agro-forest" compared to the maximum features, while there is no clear difference in performance between VV and VH. Consequently, the minimum and mean features were found to be more sensitive to mining-induced changes than the maximum. Table 1. Results of the distance (DIST) and similarity (SIM) metrics of classes "water" and "agro-forest" for the investigated time series features minimum (MIN), maximum (MAX), mean (MEAN) and for both polarization channels (VV/VH) (see Figure 6). The two metrics quantify the separation between the frequency distributions of both classes. Higher values in DIST and lower values in SIM point to a better separation offered by the feature under observation. VH polarizations. Analysis of the histograms indicate similarities in the statistics derived from the minimum and mean difference in intensity for both changed and unchanged areas. For example, in the case of changed areas, an average value of +4.6 dB (standard deviation of 2.3 dB) is observed for VV minimum and a corresponding value of +4.1 dB (standard deviation of 1.8 dB) was observed for VV mean. Corresponding average values of unchanged areas are +0.01 dB (VV minimum) and +0.32 dB (VV mean), while standard deviations are 0.88 dB (VV minimum) and 0.51 dB (VV mean). Similar observations were made for VH minimum and mean intensity plots. For changed areas an average value of +5.3 dB (standard deviation of 2.3 dB) was observed in the VH minimum and an average of +4.7 dB (standard deviation of 2.1 dB) was observed in the VH mean intensity. Corresponding values for unchanged areas are −0.01 dB (VH minimum) and +0.08 dB (VH mean), while standard deviations are 0.82 dB (VH minimum) and 0.47 dB (VH mean). This indicates that changed and unchanged areas exhibit similar characteristics in both VV and VH polarizations irrespective of the statistic (minimum or mean) being used to detect difference in intensity from one period to the other. Nonetheless, in all cases, changed areas have a higher variance and intensity difference than unchanged areas. Table 2 presents the threshold values for detecting changed areas from the statistics (minimum and mean), as determined from the intersection of the respective histograms of changed and unchanged areas/targets ( Figure 6) as well as from the Otsu algorithm. For VH polarization, the Difference images of the VV/VH minimum/mean intensities between periods "2015/2016" and "2016/2017" were used to derive thresholds for detecting changes in illegal mining areas. Results of this analysis are shown in Figure 6 and Table 2. The high-resolution Google Earth images (see Figure 3), available for a smaller portion of the study area, were analyzed to identify targets ("changed"/"unchanged"), which served as the basis for evaluating the time series features minimum and mean and determining a threshold that represents/indicates change. Figure 6 shows that changed areas generally had higher difference in intensity and variance than unchanged areas in both VV and VH polarizations. Analysis of the histograms indicate similarities in the statistics derived from the minimum and mean difference in intensity for both changed and unchanged areas. For example, in the case of changed areas, an average value of +4.6 dB (standard deviation of 2.3 dB) is observed for VV minimum and a corresponding value of +4.1 dB (standard deviation of 1.8 dB) was observed for VV mean. Corresponding average values of unchanged areas are +0.01 dB (VV minimum) and +0.32 dB (VV mean), while standard deviations are 0.88 dB (VV minimum) and 0.51 dB (VV mean). Table 2. Results of the distance (DIST) and similarity (SIM) metrics for classes "unchanged" and "changed" for time series features minimum (MIN) and mean (MEAN) in both polarization channels (VV/VH). The two metrics quantify the separation between the frequency distributions of both classes (see Figure 7). Higher DIST values and lower SIM values indicate a better separation between the classes. The last two columns show the thresholds estimated by the histogram intersection and by the Otsu Algorithm. Corresponding values for unchanged areas are −0.01 dB (VH minimum) and +0.08 dB (VH mean), while standard deviations are 0.82 dB (VH minimum) and 0.47 dB (VH mean). This indicates that changed and unchanged areas exhibit similar characteristics in both VV and VH polarizations irrespective of the statistic (minimum or mean) being used to detect difference in intensity from one period to the other. Nonetheless, in all cases, changed areas have a higher variance and intensity difference than unchanged areas. Table 2 presents the threshold values for detecting changed areas from the statistics (minimum and mean), as determined from the intersection of the respective histograms of changed and unchanged areas/targets ( Figure 6) as well as from the Otsu algorithm. For VH polarization, the histogram intersection analysis revealed a threshold value of +1.22 dB based on the mean difference in intensity and a value of +1.65 dB for the minimum difference in intensity. The corresponding threshold values identified by the Otsu algorithm were generally higher: +2.8 dB (VH minimum) and +2.30 dB (VH mean). For VV polarization similar thresholds as in the VH were found in the histogram intersection, which were +1.62 dB for the minimum difference in intensity and +1.46 dB for mean difference in intensity. The corresponding threshold values identified by the Otsu algorithm were again higher, i.e., +2.34 dB (VV minimum) and +1.99 dB (VV mean). Compared to the DIST and SIM metrics, it was found that highest values in DIST and lowest values in SIM were obtained by mean features in both polarizations. On the other hand, the minimum features showed very small differences. Again, no meaningful difference between the VH and VV polarization was found, with both exhibiting very similar values. Figure 7 shows a subset of the mean and minimum features of VV and VH polarization for periods "2015/2016" and "2016/2017", the resulting difference images, and the Google Earth Imagery. The visual comparison and interpretation indicates that mining activities can be detected with positive intensity difference values of approximately >+2 dB in both time series features. While changed areas in the upper part of the subset were clearly indicated by significant intensity differences, some areas in the lower part of the subset were not indicated.

Change Detection
Using the thresholds identified in Table 2, change maps (i.e., binary classifications of "unchanged" and "changed") were generated for all three periods (i.e., for "2015/2016"-"2016/2017", "2016/2017"-"2017/2018", and "2017/2018"-"2018/2019"). This operation is illustrated in Figure 8 for the difference in VH minimum intensities and for thresholds identified by the histogram intersection (threshold = +1.65 dB) and the Otsu algorithm (threshold = +2.78 dB). Comparing the resulting classifications, it is clear that less "changed" areas were identified by the Otsu algorithm threshold, as a higher threshold yielded a stricter prerequisite for classification of "changed" areas. Further, it is clearly visible in the imagery that in Period 1 much more changed areas were identified compared to the subsequent periods. Subsequently, accuracy assessment was carried out in order to judge, which combination of polarization, feature and thresholds finally offered the best classification accuracy based on the reference data. This assessment (Table 3) was done for class "changed" based on the differences between the VV/VH mean/minimum intensities between periods "2015/2016" and "2016/2017" (Period 1), the two thresholds identified in Table 2 and the reference data information taken from the Google Earth Imagery (Figure 3).

Change Detection
Using the thresholds identified in Table 2, change maps (i.e., binary classifications of "unchanged" and "changed") were generated for all three periods (i.e., for "2015/2016"-"2016/2017", "2016/2017"-"2017/2018", and "2017/2018"-"2018/2019"). This operation is illustrated in Figure 8 for the difference in VH minimum intensities and for thresholds identified by the histogram intersection (threshold = +1.65 dB) and the Otsu algorithm (threshold = +2.78 dB). Comparing the resulting classifications, it is clear that less "changed" areas were identified by the Otsu algorithm threshold, as a higher threshold yielded a stricter prerequisite for classification of "changed" areas. Further, it is clearly visible in the imagery that in Period 1 much more changed areas were identified compared to the subsequent periods. Subsequently, accuracy assessment was carried out in order to judge, which combination of polarization, feature and thresholds finally offered the best classification accuracy based on the reference data. This assessment (Table 3) was done for class "changed" based on the differences between the VV/VH mean/minimum intensities between periods "2015/2016" and "2016/2017" (Period 1), the two thresholds identified in Table 2 and the reference data information taken from the Google Earth Imagery (Figure 3).  Table 2) and (h, j) classification using a threshold of +2.78 dB identified via the Otsu Algorithm (see Table 2). Note: blue polylines show the channel network (see Section 2.2.2). A majority filter of 5 by 5 pixels was applied in the postclassification to classifications (e-j) (Section 2.4).  Table 2) and (h-j) classification using a threshold of +2.78 dB identified via the Otsu Algorithm (see Table 2). Note: blue polylines show the channel network (see Section 2.2.2). A majority filter of 5 by 5 pixels was applied in the post-classification to classifications (e-j) (Section 2.3.2). Table 3. Accuracy assessments of class "changed" based on the differences between the VV/VH mean/minimum intensities between sets "2015/2016" and "2016/2017" (Period 1) (see Figure 9), the thresholds identified in Table 2 and reference data generated from the Google Earth Imagery (see Figures 1 and 8). Note: Prod. Acc. = Producer Accuracy, User Acc. = User Accuracy. The assessment (Table 3) indicates that most accurate classifications (i.e., highest user accuracy (UA) and producer accuracy (PA)) were achieved in VH polarization using the minimum (UA = 72.39%, PA = 84.89%) or mean (UA = 83.55%, PA = 77.10%) features and the threshold identified by the histogram intersection. The stricter (higher) threshold values of the Otsu algorithm generally lead to higher UAs at the expense of PA (i.e., they are indicated by lower commission but higher omission). Therefore, best UA of 92.58% is found for the VH mean intensities and the threshold value of the Otsu algorithm; however, this classification also shows the second lowest PA with 46.48%, lowest PA with 37.13 % is observed for VV mean intensities and the threshold identified by the Otsu algorithm. Based on these results, the final classifications were processed using the VH polarization, the minimum feature and the threshold identified via the histogram intersection. The results are shown in Figure 9. The classifications are further displayed as density maps. These were processed by summing up changed pixels within a radius of 1 km ( Figure 10). The results indicate that land cover changes in Period 1 were frequently found along the Ofin, Oda and Pra Rivers, as well as along the Ankobra River and its tributaries in the southwestern part of the study area. The change patterns were clearly oriented and followed the course of the channel network for most locations. Generally, less land cover changes were found in Period 2. Here, the highest density of changes were found along the Ofin River and the Ankobra River south of Bekwai. Similarly, changes in Period 3 were less frequent compared to the preceding periods and highest density was found along the Ofin River and the Ankobra River south of Bekwai.

Feature
These observations are summarized in Figure 11. The total changed area reduced over time from about 102 km 2 ("2015/2016" and "2016/2017"), to 60 km 2 ("2016/2017" and "2017/2018"), to 33 km 2 ("2017/2018" and "2018/2019"). This means that the changed areas halved from first to second, and from second to third periods. Along with this, it was found that the average Shreve stream order of the changed areas reduced over time (Figure 11b), which indicates that activities moved to smaller tributaries in the course of time. This decrease was strongest from first to second period, while changes in Shreve order were lower between second and third periods. Considering the spatial patterns of changed areas (Figures 9 and 10) this was related to fewer changes observed along the Oda and Pra Rivers, as these river systems were characterized by high values of the Shreve stream order.

Susceptibility of Sentinel-1 data to Cloud Cover
The analysis revealed that a significant portion of the S-1 data were affected by atmospheric effects, typical for SAR acquisitions taken over cloud-prone regions [37]. Affected areas were characterized by a strong local decrease in the backscatter intensity of single scenes compared to the long-term mean intensity of the time-series. It is assumed that the strong lowering of the backscatter was caused by large convection cells and heavy rainfall events that are typical of the study region. However, this effect is usually more frequently observed at short wavelengths (e.g., X-Band) [37,38], while longer wavelength such as C-and L-Bands are usually less affected due to the larger penetration depth. Nevertheless, Alpers et al. [38] pointed out that high precipitation amounts of several tens of millimeters per hour can cause a significant lowering of the backscattering intensity, even at C-Band. The S-1 scenes acquired in March or October (i.e., during the rainy seasons) were more frequently affected than scenes acquired in other months ( Figure 12). The first (major) season begins in March, while the second (minor) usually starts in September [39]. Especially at the beginning of the first season and the end of the second season (i.e., the months of March/April and October/November) convective rainfalls, such as strong local rainstorms and squall lines (i.e., line of thunderstorms), are frequently observed in the area. These events are characterized by intensive and short-lasting rainfalls, which can reach precipitation amounts of more than 30 mm/h within the first half hour [39]. Friesen [40] noted that rainstorms are short lasting (1-2 hours) and stationary systems that usually affect rather small areas of approx. 20-50 km²; a size that fits well with the average extent of the affected areas. Further, the rainstorm events occur most frequently in the evening due to their convective nature and this is also the time at which the S-1 images are acquired.

Susceptibility of Sentinel-1 data to Cloud Cover
The analysis revealed that a significant portion of the S-1 data were affected by atmospheric effects, typical for SAR acquisitions taken over cloud-prone regions [37]. Affected areas were characterized by a strong local decrease in the backscatter intensity of single scenes compared to the long-term mean intensity of the time-series. It is assumed that the strong lowering of the backscatter was caused by large convection cells and heavy rainfall events that are typical of the study region. However, this effect is usually more frequently observed at short wavelengths (e.g., X-Band) [37,38], while longer wavelength such as C-and L-Bands are usually less affected due to the larger penetration depth. Nevertheless, Alpers et al. [38] pointed out that high precipitation amounts of several tens of millimeters per hour can cause a significant lowering of the backscattering intensity, even at C-Band. The S-1 scenes acquired in March or October (i.e., during the rainy seasons) were more frequently affected than scenes acquired in other months (Figure 12). The first (major) season begins in March, while the second (minor) usually starts in September [39]. Especially at the beginning of the first season and the end of the second season (i.e., the months of March/April and October/November) convective rainfalls, such as strong local rainstorms and squall lines (i.e., line of thunderstorms), are frequently observed in the area. These events are characterized by intensive and short-lasting rainfalls, which can reach precipitation amounts of more than 30 mm/h within the first half hour [39]. Friesen [40] noted that rainstorms are short lasting (1-2 h) and stationary systems that usually affect rather small areas of approx. 20-50 km 2 ; a size that fits well with the average extent of the affected areas. Further, the rainstorm events occur most frequently in the evening due to their convective nature and this is also the time at which the S-1 images are acquired. Remote Sens. 2019, 11, x FOR PEER REVIEW 20 of 26

Change Detection and Monitoring of Illegal Mining
This study employed a change detection method that was based on the generation of binary changed/unchanged masks using single thresholds, estimated via histogram analysis and the Otsu algorithm. The analyses were based on S-1 time-series features (i.e., minimum, mean and maximum backscattering intensity observed within one year), as they offered the following advantages compared to the analysis of single scenes: (i) speckle is reduced due to the temporal averaging, (ii) seasonal effects are averaged and compensated and thus do not appear as change areas in the classification result, (iii) the temporal features still preserve most important information on the backscatter characteristics and therefore on potential LULC changes. Due to the high signal stability of S-1 [32,33] it was expected that significant changes in the average backscattering intensities (approx. > +/-1 dB) between the time series features of two temporal periods can be attributed to changes in land cover.
As pointed out, mining activities in the study area usually cause a removal of vegetation and the creation of dug-outs. The analysis pointed out that LULC changes associated with mining activities were indicated by changes in both VV and VH intensities of +1.65 dB or higher (i.e., threshold estimated for VH minimum via histogram intersection), this compares well with the expected direction and magnitude of a land cover change from vegetation to water/bare soil. On the other hand, unchanged areas showed average values around +/-0 dB and standard deviations less than 0.9 dB for both polarizations and for all time series features. This observation is in line with results proposed by References [31,32]. Note in this context that the design of the approach inferred that only new sites could be detected, e.g., no information on the activity of already established sites or their abandonment was provided.
The analysis of the distance metrics (SIM and DIST), as well as the accuracy assessments conducted for Period 1 and for a smaller subset of the study area, showed that VV and VH polarization perform similar in change detection. However, the VH polarization achieved better classification accuracies, which might be attributed to the higher sensitivity of VH polarization to volume scattering processes. This sensitivity in turn leads to a better separation of water and vegetation classes in the feature space, compared to the VV polarization. Scattering from water is characterized by very low volume scattering contribution, whereas the presence of vegetation causes significantly higher volume scattering contributions and intensities at C-Band [41,42].
Among the investigated time series features, the minimum feature performed better than the mean or the maximum feature. This is reasonable, as the minimum feature is sensitive to any changes in backscattering within the period, as just the lowest intensity value observed is displayed in this

Change Detection and Monitoring of Illegal Mining
This study employed a change detection method that was based on the generation of binary changed/unchanged masks using single thresholds, estimated via histogram analysis and the Otsu algorithm. The analyses were based on S-1 time-series features (i.e., minimum, mean and maximum backscattering intensity observed within one year), as they offered the following advantages compared to the analysis of single scenes: (i) speckle is reduced due to the temporal averaging, (ii) seasonal effects are averaged and compensated and thus do not appear as change areas in the classification result, (iii) the temporal features still preserve most important information on the backscatter characteristics and therefore on potential LULC changes. Due to the high signal stability of S-1 [32,33] it was expected that significant changes in the average backscattering intensities (approx. > ±1 dB) between the time series features of two temporal periods can be attributed to changes in land cover.
As pointed out, mining activities in the study area usually cause a removal of vegetation and the creation of dug-outs. The analysis pointed out that LULC changes associated with mining activities were indicated by changes in both VV and VH intensities of +1.65 dB or higher (i.e., threshold estimated for VH minimum via histogram intersection), this compares well with the expected direction and magnitude of a land cover change from vegetation to water/bare soil. On the other hand, unchanged areas showed average values around ±0 dB and standard deviations less than 0.9 dB for both polarizations and for all time series features. This observation is in line with results proposed by References [31,32]. Note in this context that the design of the approach inferred that only new sites could be detected, e.g., no information on the activity of already established sites or their abandonment was provided.
The analysis of the distance metrics (SIM and DIST), as well as the accuracy assessments conducted for Period 1 and for a smaller subset of the study area, showed that VV and VH polarization perform similar in change detection. However, the VH polarization achieved better classification accuracies, which might be attributed to the higher sensitivity of VH polarization to volume scattering processes. This sensitivity in turn leads to a better separation of water and vegetation classes in the feature space, compared to the VV polarization. Scattering from water is characterized by very low volume scattering contribution, whereas the presence of vegetation causes significantly higher volume scattering contributions and intensities at C-Band [41,42].
Among the investigated time series features, the minimum feature performed better than the mean or the maximum feature. This is reasonable, as the minimum feature is sensitive to any changes in backscattering within the period, as just the lowest intensity value observed is displayed in this feature. Due to this fact, it was of high importance to screen and remove datasets that are affected by atmospheric conditions (i.e., clouds), as they might cause false detections. Conversely, the mean feature is dependent on the time at which the changes occur, as all intensity values of the entire period are averaged, therefore changes that occur at the end of the observation period are more difficult to detect.
It is clear that the temporal integration and the analysis of time series features has drawbacks. First, no information on the exact timing of the change is provided beside the year/period. Second, the approach is not suited to the dynamic detection of new illegal mining areas, as the approach essentially requires longer observation periods to generate the time series features. These limitations also infer that small and short-lasting changes are unlikely to be detected, also considering the spatial resolution of the S-1 data. Since this study focused on detecting inter-annual changes in rather large-scale illegal mining activities in the past, it is expected that the stated drawbacks did not have a significant effect on our results.
In this study, the thresholds for generating the binary changed/unchanged masks were estimated using a supervised (histogram intersection) and unsupervised (Otsu algorithm) approach. Both approaches showed acceptable to good classification accuracies, which were ca. 68% (Prod. Acc.) and 83% (User Acc.) for the threshold estimated via the Otsu algorithm and ca. 84% (Prod. Acc.) and 72% (User. Acc.) for the threshold estimated via the histogram intersection. Higher thresholds values were generally observed via the Otsu algorithm, which resulted in higher user but lower producer accuracies. This means that higher threshold values result in more strict and conservative estimates and have fewer false alarms at the expense of more missed detections. Higher user accuracy in turn might be a prerequisite for actual ground validation and legal actions and thus desired in an operational context.
A challenge of this study was the quantity and extent of the reference data. As pointed out, the results were validated for a smaller region using high-resolution Google Earth Imagery. Mining activities can be detected on these imagery with a high level of confidence. However, due to persistent cloud cover, only few images were available that (i) covered a sufficiently large area affected by mining activities and (ii) the period of acquisition corresponded to the S-1 observation periods making them suitable for validation. Consequently, only few images provided a useful and reliable reference that was suited to assess the classification accuracies provided by the unsupervised and supervised binary classifications. In this regard, upcoming studies should compile more suited high-resolution optical space-borne imagery, also from commercial sources, e.g., the archive of the Planet satellites (planet.com). Having more high-quality reference data will allow for a more comprehensive classification and referencing, e.g., also for the most recent observation periods. On the other hand, future studies can embark on an extensive field campaign (including historical land use interviews) to generate more reference data for validation. Considering the increased observation frequency of optical space-borne systems (e.g., due to the operation of Sentinel-2), it might also be possible to fuse SAR and optical data in the classification process in order to enhance the classification accuracies, to increase the sensitivity of the detection and the spatial resolution of the classification result.
The presented approach does not require a comprehensive processing and therefore it can be repeated by using the Ground Range Detected (GRD) products of Sentinel-1. These products require less disk space, are already partially radiometrically corrected and the processing of sigma or gamma nought requires less processing steps, compared to the SLC products. Further, it is suggested to run the approach using cloud-processing facilities such as the Google Earth Engine (GEE) (earthengine.google.com) [43]. GEE already offers the necessary infrastructure to assess the S-1 archive (i.e., the GRD products are available) and to process time series features from a stack of S-1 images. This in turn will save the download and processing of large amounts of S-1 data and will provide analyses ready data. However, the issue that S-1 data are affected by atmospheric conditions (i.e., clouds) must be addressed, as patches of low backscattering will manifest in the minimum time series feature. Therefore, other metrics, e.g., such as the analyses of lower quantile values instead of the minimum, should be investigated as a cloud-based processing strategy is not suited to the manual selection and removal of individual scenes.

Trend in the Extent of Illegal Mining Areas between 2015 and 2019
Our results showed a decreasing trend in the extent of illegal mining in the study area from 102 km 2 in Period 1 (2015/2016-2016/2017) to 60 km 2 and 33 km 2 in periods 2 (2016/2017-2017/2018) and 3 (2017/2018-2018/2019), respectively. Despite its widespread practice and resultant degradation and pollution, very few studies have used remote sensing data to quantify changes (over a period) in the extent of illegal mining activities in the country. This is principally attributable to excessive cloud cover in the southern part of Ghana where illegal mining mostly takes place, and the difficulty in obtaining usable optical images. Owusu-Nimo et al. [44] mapped and analyzed the spatial distribution patterns of artisanal small-scale gold mining in the Western region of Ghana, but did not use satellite imagery nor quantified the extent/area of the activity in the region. They generally mapped galamsey sites as points using field-based surveys and found a total of 911 galamsey sightings comprising of 7470 individual operations in 312 communities across 11 municipal and district assemblies in Western region. Out of 11 districts, they found three as the galamsey hotspot districts: Tarkwa Nsuaem (294 sightings), Amenfi East (223 sightings) and Prestea Huni Valley Districts (156 sightings). Two out of these three hotspot districts fall within our study area. Awotwi et al. [45] conducted a general land use and land cover mapping in the Pra river basin using a time-series of Landsat images acquired between 1986 and 2016. They noted an increase of 10.62% in the extent of mining areas between 2004 and 2008, and 304% between 2004 and 2016. Other studies also inferred increases in the extent of mining areas based on general land use and land cover analysis based on satellite imagery [20,46].
Snapir et al. [21] is the only study we have found that had a specific focus on assessing the extent of illegal mining sites between 2011 and 2015 using optical satellite imagery. For a subset of their study area (change area), they found that the extent of illegal mining areas tripled between 2011 and 2015, increasing from 109.07 km 2 in 2011 to 232.82 km 2 and 366.96 km 2 in 2013 and 2015, respectively. Although the area and period investigated are not the same as that of this study's, a note-worthy difference in their results compared to ours is the trend in the extent of illegal mining areas, i.e., increasing trend in Snapir et al. [21] and decreasing trend in this study. Two possible reasons for this variance can be discussed. First, as alluded by Snapir et al. [21], their approach did not differentiate between active and abandoned pits or mining areas. Thus, the area detected in, for example 2013, could be a summation of abandoned pits from 2011 and active/new pits that opened in 2013. Similarly, active/new pits in 2015 would add up to those of 2011 and 2013 to indicate a much larger extent compared to previous years. The approach adopted in this study, however, detects mining-induced changes from one period to the other. Thus, the detected changes at Period 2, relative to Period 1, can be considered as new illegal mining areas that opened in the period under analysis. In this sense, when areas mapped in, for example, Period 1 are added to the area of detected change in Period 2, there would be an overall increase in the extent of illegal mining areas at Period 2. When analyzed this way, our results will indicate an increasing trend in the extent of illegal mining areas as 102 km 2 , 162 km 2 and 195 km 2 for Periods 1, 2 and 3, respectively.
The second possible reason for the observed difference in trends, specifically the decreasing trend in our results (with respect to new illegal mining areas), is the government's efforts to sanitize the small-scale mining sector and reduce illegal mining activities in the country. Although governmental effort to curb the menace started in 2013 with the establishment of an inter-ministerial task force, efforts from 2017 seem to have been more effective. These efforts include: (i) a six-month ban on all small-scale mining activities (even including licensed operators) and seizure of equipment, (ii) suspension of the issuance of licenses for small-scale mining, (iii) setting up of a Multilateral Mining Integrated Project (MMIP) to enact more stringent mining regulations and enforcement of mining laws, (iv) designation of 14 courts by the then Chief Justice to handle galamsey cases and (v) setting up of a media coalition against galamsey to support government's efforts. In addition, an Afrobarometer survey of 2017 revealed overwhelming support of Ghanaians for government's efforts, with 74% indicating that no citizen should be permitted to engage in illegal mining for any reason [47]. These factors may have accounted for the observed decreasing trend in the extent of new illegal mining areas between 2016/2017-2017/2018 and 2017/2018-2018/2019. Our results indicate that government's efforts between 2017 and 2019 resulted in lesser areas being degraded post 2017, compared to pre-2017 periods. This notwithstanding, there is still evidence of illegal mining activities in the study area, which must be continually monitored to assess efforts at curbing the menace.

Conclusions
In this study, annual S-1 (Sentinel-1) time-series data from 2015 to 2019, together with reference data from Google Earth and existing land cover data, were analyzed using change detection techniques to map and monitor illegal mining activities on an annual basis in South-Western Ghana. The methodology was premised on the fact that mining-induced land cover changes can be detected based on changes in the scattering characteristics of time-series images acquired with the same acquisition geometry. The sensitivity of three time-series features (i.e., minimum, mean and maximum backscattering intensity observed within one year) to mining-induced changes (e.g., clearing of vegetation) was first conducted, after which analysis was performed to determine a threshold suitable to identify/classify changed and unchanged areas.
The use of time-series features proved useful for detecting inter-annual mining-induced changes along river bodies. Compared to the mean and maximum features, the minimum feature (in both VH and VV polarizations) performed better and was found to be more sensitive to changes in backscattering within the period investigated, including the lowest intensity values observed. A backscatter threshold value of +1.65 dB was found suitable to detect illegal mining activities in the study area. Application of the threshold revealed illegal mining area extents of 102 km 2 , 60 km 2  Despite the advantages of Synthetic Aperture Radar data in monitoring phenomena in cloud-prone areas, our analysis revealed reasonable susceptibility of S-1 data to atmospheric conditions due to high intensity rainfall events. Scenes acquired in March and October, representing the beginning and ending of the rainy season, were the worst affected. A total of 29 out of 117 (ca. 25%) scenes could not be used in our analysis due to a strong local decrease in backscatter in the order >−10 dB compared to the long-term mean intensity and the undisturbed surroundings. This could have negatively affected our results considering the methodology employed. Further investigation of this finding in other geographies and climatic regions is required to ascertain the suitability of S-1 data for time-series dependent analysis such as the one employed in this study. Future studies can also benefit from an integration of Sentinel-2 data into the analysis.
The quantity and quality of the reference data was found to be a major challenge of this study. Google Earth Imagery permits a detailed view of changes in illegally mined areas. However, lack of sufficient images (spatially and temporarily) over the area analyzed due to excessive cloud cover limited the amount of reference data generated for this work. Future studies should consider embarking on an extensive field data collection to map changed and unchanged areas to serve as reference data.