Simulated Annealing Applied to Satellite Data for Agricultural Monitoring

Marian Keane
Version 1
Published in
7 min readFeb 5, 2021
Photo by NASA on Unsplash

Earth observation satellites have been in orbit since the early 1970s. For decades now, public and private organisations around the world have used satellite data to monitor agricultural land in order to:

- Forecast crop production and hence commodity prices,

- Monitor risks to food supply and food security,

- Support farmers via monetary payments for approved farming practices.

Publicly-owned satellites provide freely available imagery down to 10-metre spatial resolution. More granular data can be purchased from commercial satellite owners, who can provide imagery down to 0.25 metres spatial resolution.

Countries within the European Union are preparing for the use of satellite data (“Checks by Monitoring”) to monitor eligibility criteria for income support to farmers under the Common Agricultural Policy (CAP). Imagery from the Copernicus Sentinel satellites is freely available, and member states can choose to use this data or “other data with at least equivalent value”. Machine-learning algorithms are developed and shared across member states — these algorithms use biophysical indicators derived from satellite data to identify agricultural events (ploughing, harvesting etc), classify crop type or land use, or monitor land cover change over time.

One field, one crop

The unit of land for EU CAP income support is an Agricultural parcel, defined as “a continuous area of land, declared by one farmer, which includes no more than one crop group”. If the boundary of the parcel is reliably known in advance, then machine-learning algorithms can confidently be applied to the satellite imagery, cropped to the boundaries of the parcel.

However, what if the parcel boundary is incorrect? An area of land that in the past contained only a single crop may in the current year contain two or more crops. If this is not detected, the algorithms that identify the most likely crop within the parcel boundary will surely fail.

The technical guidance to the rollout of Checks by Monitoring across the EU require that the spatial integrity of the parcel be checked, as a pre-cursor to the application of other algorithms. No specific machine-learning algorithms are prescribed, but the problem-type is that of spatial clustering — using individual pixel-level data to split an image into two or more clusters.

The plots below show NDVI (Normalised Difference Vegetation Index) for an individual parcel on sample image capture dates. Higher values of NDVI indicate the presence of live green vegetation on the ground. The April images strongly suggest the presence of two distinct crops, however the 17th April image is missing some pixel values due to cloud cover on that date. By 21st July no noticeable NDVI differences are apparent across the parcel — illustrating that any spatial clustering approach may work better on some dates than on others.

Figure 1: NDVI by pixel by date for Parcel_D

Spatial Clustering

A potential approach to identifying the locations of crop #1 and crop #2 within a 2-crop parcel is to iteratively nominate pairs of pixels as the centroids of cluster 1 and cluster 2 (where each cluster represents a distinct crop). The optimal solution on a given date is the pair of pixels that maximises the difference between the mean NDVI between the two clusters.

In the sample parcel (“Parcel_D”) above, we have up to 350 pixels in each image, each pixel representing an area of 10m². There are therefore 61,075 pairs of pixels to review (the number of ways we can choose 2 from 350).

For each such pair (Centroid1, Centroid2), the following calculations are run:

- Calculate (straight-line) distance of each other pixel to Centroid1 and Centroid2

- If a pixel is closer to Centroid1, assign it to Cluster1, otherwise to Cluster2

- Calculate the mean NDVI across each of Cluster1 and Cluster2

- Record the difference between the mean NDVIs.

Following an exhaustive search of the state-space, the pair of centroids that maximises the difference in mean NDVI is retained for each image date. These optimal centroids across a growing season are illustrated in the diagram by the blue and red points.

Figure 2: Potential centroids of Cluster 1, Cluster 2

The final weighted classification of pixels to Cluster1 or Cluster2 across the entire season is plotted below, using blue and red to indicate high probability of belonging to Cluster1 or Cluster2, with pink/purple shades indicating uncertain areas.

The results show a reasonably clear delineation of Parcel_D into two distinct areas, each representing a single crop. Now that these sub-parcels have been identified, the crop classification and other machine-learning algorithms can be safely applied to each of the sub-parcels.

Figure 3: Full-search segmentation of Parcel_D

Problem

The problem with the above spatial clustering technique is the compute time to check all pairs of pixels. The complexity of the problem is O(n²) which is simply not scalable. Parcel_D illustrated above has an area of 3.5 hectares. A 20-hectare parcel, with 2,000 pixels, would require 1,999,000 pairs of pixels (i.e. 2000 choose 2) to be checked for each image capture date. Running such calculations across all parcels to be monitored by a member state would grind an agricultural monitoring system to a halt.

Solution

One solution to the problem is to use optimisation techniques to more efficiently find the best centroids without needing to exhaustively search through all permutations of pixel-pairs.

One such optimisation technique is Simulated Annealing. Simulated annealing algorithms are based on mimicking the thermodynamic forces observed when a metal or other material is heated and then cooled. The algorithm initially runs ‘hot’, jumping around the solution space to find global maxima or minima, and over time ‘cools down’ to explores solutions more local to the current solution.

The purpose of the “temperature” feature is to avoid getting ‘stuck’ in a sub-section of the state space, since “transitions out of a local optimum are always possible at nonzero temperature” (Kirkpatrick, S.; Gelatt Jr, C. D.; Vecchi, M. P. (1983). “Optimization by Simulated Annealing”).

Figure 4: Example of a state space with a global and local maximum

The algorithm starts from an initial state s0 and continues until a maximum of kmax steps have been taken. At each iteration, the algorithm considers some neighbouring state of the current state, and probabilistically decides between moving to the neighbouring state or staying in the current state. If the temperature is high, the probability of moving to the neighbouring state is higher (even if the neighbouring state produces a worse solution). As the temperature reduces, the probability of moving to the neighbouring state reduces.

A Wikipedia entry for simulated annealing provides pseudo-code for the method:

Figure 5: Pseudo-code for simulated annealing

Each of the steps in the pseudo-code needs to be tailored to the problem at hand. Here is how the method was applied to the 2-crop problem:

Figure 6: Application of simulated annealing pseudo-code to 2-crop problem

Results

If the method works well, the identification of sub-parcels using Simulated Annealing should be reasonably close to the “full-search” of the state space, and the algorithm should run much more quickly.

Dealing first with the run-time, the simulated annealing approach was found to run at multiples of the speed of the full-search approach, particularly on larger parcels:

Figure 7: Run-time for Simulated Annealing vs Full-search

In terms of the quality or accuracy of the final solution, the difference in mean NDVI between Cluster1 and Cluster2 under the Simulated Annealing was at times lower than the maximum found using the full search method. On investigation, it was found that the images on which the Simulated Annealing approach had performed poorly tended to have issues with missing pixels (due to cloud cover). These images tended to result in the full-search method finding a “niche” solution that was an outlier relative to neighbouring pixel pairs. Increasing the value of kmax to allow the Simulated Annealing to run for longer did not materially improve the results for these cloud-affected images.

Figure 8: Full-search segmentation of Parcel_D
Figure 9: Simulated Annealing segmentation of Parcel_D

In practical terms, the inability of the Simulated Annealing approach to always return the absolute optimal result was felt to be less important than producing a reasonably good result in optimal time. Comparing the full-search and Simulated Annealing images in Figures 8 and 9 above, it can be seen that the simulated annealing returns a small area where the segmentation is less certain (purple shading) close to the boundary between the two crops, but these purple pixels are still identified as more similar to the “blue” crop than the “red”.

Conclusion

Spatial clustering algorithms, if not carefully designed, can be very slow to run. In the context of a national monitoring of agricultural activity, an algorithm with complexity of O(n²) is not practical. Simulated Annealing can deliver excellent results in a small fraction of the run-time. The algorithm requires some effort to tailor the generic method to the problem at hand, but the savings in run-time mean that the effort is well worthwhile.

--

--