Bridging the last mile in crop-type mapping

Bridging the last mile in crop-type mapping

📆 Project Period	February, 2022
📍 GitHub	gitlab.esa.int

Project Summary

This is an aggregation of different crop types in several regions of the world. The variance in the data should drive the development of new methods and different experiments. The merging of data represents only one part of the project. For each field polygon, a time series of mean values and other statistical properties was created. For example, a field (polygon) is represented by a Sentinel-2 time series of means per field. The labels for the crop types come from different sources and have been summarised.

Notebook

Example for crop type data

Crop type labels and id:

Bavaria: Potatos: 600, 601, 602 Winter barley: 476, 131 Winter rapeseed: 311, 489 Winter wheat: 115 Corn: 400, 177, 171, 411, 410 Sugar beet: 603

The labels for the following data are documented in:

Kenya: https://mlhub.earth/data/ref_african_crops_kenya_01
Central Asia: https://www.nature.com/articles/s41597-020-00591-2?proof=t
France: https://breizhcrops.org/

#plot shapefile for durnast in bavaria
durnast_shp.plot()

png

png

png

Crop type labels

durnast_ts[durnast_ts.id == '58cb970d59ab4e618787ff1aadd85599'].crop_type.unique()

References

[1] Breizhcrops, https://breizhcrops.org

[2] Radiant Earth Foundation, https://www.radiant.earth/

[3] Remelgado, R., Zaitov, S., Kenjabaev, S. et al. A crop type dataset for consistent land cover classification in Central Asia. Sci Data 7, 250 (2020). https://doi.org/10.1038/s41597-020-00591-2

[4] Chair of Plant Nutrition, TUM, https://www.pe.wzw.tum.de/

# Demo for visualisation of crop type and yield data
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LinearRegression

#3D stuff
from IPython.core.display import display, HTML
import json

import geopandas as gpd
import geojson
import pandas as pd

from functools import partial
from shapely.geometry import Point # Point class
from shapely.geometry import shape # shape() is a function to convert geo objects through the interface

#!pip install matplotlib
import matplotlib.pyplot as plt

#load shapefiles
centralasia_shp = gpd.read_file(
    "../data/cropdata/CentralAsia/CAWa_CropType_samples.shp")

durnast_shp = gpd.read_file(
    "../data/cropdata/Bavaria/bavaria_shs/labels_new.shp")

kenya_shp = gpd.read_file(
    '../data/cropdata/Kenya/kenya_shs/labels_new.shp')

# load sentinel-2 timeseries L2A for shapefiles
# The timeseries are from GEE and Sentinel Hub
# Sentinel hub provides a FIS API and a new SHS API. Both data sets are included
# the timeseries for both levels L1C and L2A are available
# <https://sentinelhub-py.readthedocs.io/en/latest/examples/fis_request.html>
# all data from Sentinel hub is cloud free

durnast_ts = pd.read_excel(
    "../data/cropdata/Bavaria/bavaria_shs/timeseries_L2A.xlsx")

kenya_ts = pd.read_excel(
    '../data/cropdata/Kenya/kenya_shs/timeseries_L2A.xlsx')

bavaria_train = pd.read_excel(
    "../data/cropdata/Bavaria/sentinel-2/Training_bavaria.xlsx")
bavaria_test = pd.read_excel(
    "../data/cropdata/Bavaria/sentinel-2/Test_bavaria.xlsx")
bavaria_reordered = pd.read_excel(
    '../data/cropdata/Bavaria/sentinel-2/data2016-2018.xlsx', index_col=0)
bavaria_test_reordered = pd.read_excel(
    '../data/cropdata/Bavaria/sentinel-2/TestData.xlsx', index_col=0)

#clean up
durnast_ts.drop(['Unnamed: 0'], axis=1, inplace=True)
bavaria_test.drop(['Unnamed: 0'], axis=1, inplace=True)
bavaria_train.drop(['Unnamed: 0'], axis=1, inplace=True)
bavaria_train.Date = bavaria_train.Date.dt.strftime('%m-%d')

# plot field and corresponding timeseries L2A of Sentinel-2
# every field has a unique id, the labels are included in the shapefile and corresponding timeseries

fig, (ax1, ax2) = plt.subplots(1, 2,figsize=(12, 5))
durnast_ts[durnast_ts.id == '58cb970d59ab4e618787ff1aadd85599'][['bands_B03_mean','bands_B12_mean']].plot(ax=ax2)
durnast_shp[durnast_shp.id == '58cb970d59ab4e618787ff1aadd85599'].plot(ax=ax1)

#plot corresponding NDVI timeseries L2A for a field
fig, ax = plt.subplots(figsize=(15, 8))

series = durnast_ts[durnast_ts.id == '58cb970d59ab4e618787ff1aadd85599']

series.plot(
    ax=ax,
    x='interval_from',
    y='bands_NDVI_mean')

ax.fill_between(
    series.interval_from.values,
    series['bands_NDVI_mean'] - series['bands_NDVI_stDev'],
    series['bands_NDVI_mean'] + series['bands_NDVI_stDev'],
    alpha=0.3
)