Segment buildings in African cities from aerial imagery and advance Responsible AI ideas for disaster risk management

Post by Dave Luo, Grace Doherty, Nicholas Jones, and Vivien Deparday, GFDRR Labs/World Bank

Takeaways

The Global Facility for Disaster Reduction and Recovery (GFDRR) is partnering with Azavea and DrivenData to introduce a new dataset and machine learning (ML) competition ($15,000 in total prizes) to improve mapping for resilient urban planning. Better ML-supported mapping for disaster risk management means addressing barriers to applying ML in African urban environments and adopting best practices in geospatial data preparation to enable easier ML usage. The competition dataset — over 400 square kilometers of high-resolution drone imagery and 790K building footprints — is sourced from locally validated, open source community mapping efforts from 10+ urban areas across Africa. Prize-winning solutions will be published as open-source tools for continued ML development and benchmarking.

The Open Cities AI Challenge has two participation tracks:

1. $12,000 in prizes for best open-source semantic segmentation of building footprints from drone imagery that can generalize across a diverse range of African urban environments, spatial resolutions, and imaging conditions.
2. $3,000 in prizes for thoughtful explorations of Responsible AI development and application for disaster risk management. How might we improve the creation and use of ML systems to mitigate biases, promote fair and ethical use, inform decision-making with clarity, and make safeguards to protect users and end-beneficiaries?

The competition is ongoing and ends March 16th, 2020. Join today!

Open Data for Resilient Urban Planning

Cities around the world are growing rapidly, especially in Africa — by 2030, half of Sub-Saharan Africa’s population will live in urban areas. As urban populations grow, their exposure to flooding, erosion, earthquakes, coastal storms, and other hazards becomes a complex challenge for urban planning.

Understanding how assets and people are vulnerable to these risks requires detailed, up-to-date geographic data of the built environment. For example, a building’s particular location, shape, and construction style can tell us whether it will be more exposed to earthquake or wind damage than nearby buildings. Roads, buildings, and critical infrastructure need to be mapped frequently, accurately, and in detail if we are to understand and manage risk effectively. But in countries with less developed data infrastructure, traditional urban data collection methods can’t keep up with increasing density and sprawl.

Thankfully, collaborative and open data collection practices are reshaping the way we map cities. Today, local mapping communities are improving maps for some of the world’s most vulnerable neighborhoods — bringing highly accurate and detailed geographic data up-to-date and to scale. GFDRR at the World Bank supports programs like Open Cities Africa and Dar Ramani Huria to map buildings, roads, drainage networks and more in over a dozen African cities, and Zanzibar Mapping Initiative was the world’s largest aerial mapping exercise using consumer drones and local mappers to produce open spatial data for conservation and development in the archipelago.

Data collected in these community mapping programs are used to design tools and products that support government decision-making. Digitized maps are published to OpenStreetMap and aerial imagery to OpenAerialMap where they serve as data public goods that can be used and improved by all. The open source philosophy behind the movement and an emphasis on local skill-building has fostered local networks of talent in digital cartography, robotics, software development, and data science.

Potential of Machine Learning for Mapping

Advances in ML for visual tasks could further improve mapping quality, speed, and cost. Recent examples of ML applications for mapping include Facebook’s AI-assisted mapping tool for OpenStreetMap and Microsoft’s country-scale automated building footprint extraction (in USA, Canada, Tanzania and Uganda). Competitions like SpaceNet and xView2 advance ML practices for automated mapping of roads, buildings, and building damage assessment after disasters.

Obstacles, however, stand in the way of effectively applying current ML mapping solutions to the African disaster risk management context. Africa’s urban environments differ significantly in make-up and appearance from European, American, or Asian cities which have more abundant data that ML models are often trained on.

Buildings that are more densely situated and diverse in shape, construction style, and size may be less recognizable to ML models that saw few or no such examples in their training.

Imagery is collected by commercial drones at much higher resolution under diverse environmental conditions, requiring adaptation of models usually trained on lower-resolution, more consistently collected and preprocessed satellite imagery.

Crowdsourced and community-driven data labeling may differ greatly in what base imagery layers are used, workflow, data schema, and quality control, requiring models that are robust to more label noise.

Geospatial data comes in a diversity of file formats, sizes, and schemas that create high adoption and knowledge barriers that hamper their use in machine learning.

There is now a growing abundance of locally-validated open map data and high resolution drone imagery in diverse built environments. How might we best address these obstacles and enhance the state of practice in machine learning to support mapping for urban development and risk reduction for Africa’s cities?

Introducing the Open Cities AI Challenge

Dataset

Working with partners Azavea and DrivenData, the Labs team at GFDRR combined the excellent work of many participatory mapping communities across Africa, applied best practices in cloud-native geospatial data processing (i.e. using Cloud-Optimized GeoTIFFs [COG] and SpatioTemporal Asset Catalogs [STAC]), and standardized wherever possible to make data more readily usable for machine learning. The result is a novel, extensive, open dataset of over 790K building footprints and 400 square kilometers of drone imagery representing 10 diverse African urban areas in ML-ready form.

Using COG and STAC for geospatial data provides us with bandwidth-efficient, rapid, and query-able access to our imagery and labels in a standardized format. Ease of access to files and indexing of data catalogs is particularly important for geospatial data which can quickly grow to 100s of gigabytes. It also enables us to tap into the growing ecosystem of COG and STAC tools, like STAC Browser to rapidly visualize and access any training data asset in a web browser, despite individual image files being up to several GBs and the entire dataset totaling over 70 GBs in size:

PySTAC, a new Python library by Azavea, enables STAC users to load, traverse, access, and manipulate data within catalogs programmatically. For example, reading a STAC catalog:

train1_cat = Catalog.from_file('https://drivendata-competition-building-segmentation.s3-us-west-1.amazonaws.com/train_tier_1/catalog.json') 
train1_cat.describe()

* <Catalog id=train_tier_1>
    * <Collection id=acc>
      * <Item id=665946>
      * <LabelItem id=665946-labels>
      * <Item id=a42435>
      * <LabelItem id=a42435-labels>
      * <Item id=ca041a>
      * <LabelItem id=ca041a-labels>
      * <Item id=d41d81>
      * <LabelItem id=d41d81-labels>
    * <Collection id=mon>
      * <Item id=401175>
      ...

Inspecting an item’s metadata:

one_item = train1_cat.get_child(id='acc').get_item(id='ca041a')
one_item.to_dict()

{
  "assets": {
    "image": {
      "href": "https://drivendata-competition-building-segmentation.s3-us-west-1.amazonaws.com/train_tier_1/acc/ca041a/ca041a.tif",
      "title": "GeoTIFF",
      "type": "image/tiff; application=geotiff; profile=cloud-optimized"
    }
  },
  "bbox": [
    -0.22707525357332697,
    5.585527399115482,
    -0.20581415249279408,
    5.610742610987594
  ],
  "collection": "acc",
  "geometry": {
    "coordinates": [
      [
        [
          -0.2260939759101167,
          5.607821019807083
        ],
        ...
        [
          -0.2260939759101167,
          5.607821019807083
        ]
      ]
    ],
    "type": "Polygon"
  },
  "id": "ca041a",
  "links": [
    {
      "href": "../collection.json",
      "rel": "collection",
      "type": "application/json"
    },
    {
      "href": "https://drivendata-competition-building-segmentation.s3-us-west-1.amazonaws.com/train_tier_1/acc/ca041a/ca041a.json",
      "rel": "self",
      "type": "application/json"
    },
    {
      "href": "../../catalog.json",
      "rel": "root",
      "type": "application/json"
    },
    {
      "href": "../collection.json",
      "rel": "parent",
      "type": "application/json"
    }
  ],
  "properties": {
    "area": "acc",
    "datetime": "2018-11-12 00:00:00Z",
    "license": "CC BY 4.0"
  },
  "stac_version": "0.8.1",
  "type": "Feature"
}

Learn more about the dataset and STAC resources.

Competition

Accompanying the dataset is a competitive machine learning challenge with $15,000 in total prizes to encourage ML experts globally to develop more accurate, relevant, and readily usable open-source solutions to support mapping in African cities. There are 2 participation tracks:

Semantic Segmentation track: $12,000 in prizes for the best open-source semantic segmentation models to map building footprints from aerial imagery.

The machine learning objective is to segment (classify) every pixel in every image as building or no-building with model performance being evaluated with the Intersection-over-Union metric (aka Jaccard Index):

Semantic segmentation is useful for mapping because its pixel-level outputs are relatively easy to visually interpret, verify, and use as-is (e.g. in the calculation of built-up surface area) or as inputs to downstream steps (e.g. first segment buildings and then classify attributes about each segmented building like its construction status or roof material).

Segmentation track participants must also submit at least once to the Responsible AI track to qualify for $12,000 in segmentation track prizes.

Responsible AI track: $3,000 in prizes will be awarded for best ideas applying an ethical lens to the design and use of ML systems for disaster risk management.

ML can improve data applications in disaster risk management, especially when coupled with computer vision and geospatial technologies, by providing more accurate, faster, or lower-cost approaches to assessing risk. At the same time, we urgently need to develop a better understanding of the potential for negative or unintended consequences of their use. With growing attention given to questions of appropriate and ethical ML use for facial recognition, criminal justice, healthcare, and other domains, we have an immediate responsibility to elevate these questions for disaster risk.

Examples of potential harm that ML technologies present in this space include, but are not limited to:

• Perpetuating and aggravating societal inequalities through the presence of biases throughout the machine learning development pipeline.
• Aggravating privacy and security concerns in Fragility, Conflict and Violence settings through combination of previously distinct datasets.
• Limiting opportunities for public participation in disaster risk management due to increased complexity of data products.
• Reducing the role of expert judgement in data and modeling tasks and in turn increasing probability of error or misuse.
• Inadequately communicating methods, results, or degrees of uncertainty, which increases the chance of misuse.

ML practitioners and data scientists are uniquely positioned to examine and influence the ethical implications of our work. We ask challenge participants to consider the applied ethical issues that arise in designing and using ML systems for disaster risk management. How might we improve the creation and application of ML to mitigate biases, promote fair and ethical use, inform decision-making with clarity, and make safeguards to protect users and end-beneficiaries?

This track’s submission format is flexible: participants can submit Jupyter notebooks, slides, blogs, essays, demos, product mockups, speculative fiction, art work, synthesis of research papers or original research, or whatever other format best suits you. Submissions will be evaluated by a panel of judges on thoughtfulness, relevance, innovation, and clarity.

What Comes Next

This challenge will produce new public goods that advance our state of practice in applying ML for understanding risk in urban Africa; this includes new ML performance benchmarks for building segmentation from aerial imagery in relevant geographies, top-performing solutions for mapping in African cities, and in-depth explorations of how we responsibly create and deploy AI systems for disaster risk management.

Prize-winning solutions will be published as open-source tools and knowledge and the challenge dataset will remain an open data resource for continued ML development and benchmarking. GFDRR will use lessons learned to inform policies and procurement strategies for using ML for urban mapping and planning.

Join the Challenge!

The competition is currently running until March 16, 2020. With one month to go, there is plenty of time to explore the data and participate in either tracks but don’t delay, join today at:

drivendata.org/competitions/60/building-segmentation-disaster-resilience

The Open Cities AI Challenge was originally published in TDS Archive on Medium, where people are continuing the conversation by highlighting and responding to this story.

by daveluo (on github and elsewhere)

In this post and the accompanying Google Colab notebook, we’ll learn all the code and concepts comprising a complete workflow to automatically detect and delineate building footprints (instance segmentation) from drone imagery with cutting edge deep learning models.

All you’ll need is a Google account, an internet connection, and a couple of hours to learn how to make the data & model that learns to make something like this:

In modular steps, we’ll learn to…

Preprocess image geoTIFFs and manually labeled data geoJSON files into training data for deep learning:

Create a U-net segmentation model to predict what pixels in an image represent buildings (and building-related features):

Test our model’s performance on unseen imagery with GPU or CPU:

Post-process raw model outputs into geo-registered building shapes evaluated against ground truth:

And along the way, we’ll get familiar with great geospatial data & deep learning tools/resources like:

• Geopandas: “an open source project to make working with geospatial data in python easier. GeoPandas extends the datatypes used by pandas to allow spatial operations on geometric types.”
• Rasterio: “reads and writes geospatial raster datasets”
• Supermercado: “supercharger for Mercantile” (spherical mercator tile and coordinate utilities)
• Rio-tiler: “Rasterio plugin to read mercator tiles from Cloud Optimized GeoTIFF dataset”
• Solaris: “Geospatial Machine Learning Analysis Toolkit” by Cosmiq Works
• Cloud-Optimized GeoTIFFs (COG): “An imagery format for cloud-native geospatial processing”
• Spatio-Temporal Asset Catalogs (STAC): “Enabling online search and discovery of geospatial assets”
• OpenAerialMap: “The open collection of aerial imagery”
• Fast.ai for geospatial deep learning: “The fastai library simplifies training fast and accurate neural nets using modern best practices” built on the PyTorch framework.

How to get the most out of this tutorial:

The Colab notebook is our main learning resource — working interactively within the GPU-accelerated Colab notebook environment is highly recommended!

Code is organized into modular sections, set up for installation/import of all required dependencies, and executable on either CPU or GPU runtimes (depending on the section). Links to load files generated at each step are also included so you can pick up and start from any section. Inline# comments (& references for further reading) are provided within code cells to explain steps or nuances in more detail as needed. Executing all code cells end-to-end takes <1 hour on GPU.

This Medium post serves as a high-level conceptual walkthrough and maps directly to sections within the Colab notebook. This post works best as a quick overview with handy bookmarks to Colab (see [Colab section link] under each section heading) or viewed side-by-side with the Colab notebook as a code & concept companion set.

This tutorial assumes you have a working knowledge of Python, data analysis with Pandas, making training/validation/test sets for machine learning, and a beginner practitioner’s grasp of deep learning concepts. Or the motivation to gain what knowledge you’re missing by following the ample references linked throughout this post and notebook.

With that as mental prep, let’s do some geospatial deep learning!

Open tutorial notebook and create your own working Colab copy with File > Save a Copy in Drive (recommended option):

Or preview notebook non-interactively (not recommended for optimal learning but good for a quick look):

Notebook on nbviewer

Pre-Processing

Preview and load imagery and labels

[Colab section link]

For this tutorial, we’ll use the Tanzania Open AI Challenge dataset of 7-cm resolution drone imagery and building footprint labels over Unguja Island, Zanzibar, Tanzania. Much thanks to the following organizations for producing, openly licensing, and making this invaluable dataset accessible:

• CC-BY-4.0 licensed by Commission for Lands (COLA) — Revolutionary Government of Zanzibar (RGoZ)
• Labeled data produced & processed by State University of Zanzibar (SUZA), World Bank OpenDRI, WeRobotics
• Drone imagery created by Zanzibar Mapping Initiative and hosted on OpenAerialMap:

For simplicity of demonstration, we’ll create training and validation data from a single drone image (in cloud-optimized geoTIFF format) and its accompanying ground-truth labels of manually traced building outlines (in GeoJSON format).

We’ll work with imagery and labels from image grid znz001 which covers the northern tip of Zanzibar’s main island of Unguja. Here is a browsable preview of znz001's drone imagery with its accompanying building outline labels, indexed according to the Spatio-Temporal Asset Catalog (STAC) label extension and visualized in an instance of STAC browser:

After previewing the labeled data and imagery in the browser, we’ll import our geo-processing tools, copy the direct download URLs from the Assets tab of the browser, and test loading them in our notebook.

Draw train and validation areas of interest (AOIs) with geojson.io

[link to Colab section]

Since we are working with a single image, we need to delineate what sub-areas of the image and labels should be used as training versus validation data for model training.

Using geojson.io, we’ll draw ourtrn and val Areas of Interest (AOI) polygons in geojson format and add dataset:trn or dataset:val to the respective polygon properties.

The finished polygons will look something like this in geojson.io:

And here is the exact GeoJSON file I created displayed in github gist:

For demonstration of later steps, I intentionally drew a more complex shape for each AOI but we could have simply drawn adjacent rectangles instead.

Or in more complex cases, we could choose to draw AOIs of smaller sub-areas that don’t encompass the entire image — for instance, if we want to create training data for specific types of environments like dense urban areas or sparsely populated rural areas only or we want to avoid using poorly labeled areas in our training data.

Drawing the AOIs as geoJSON polygons in this way gives us the flexibility to choose exactly what and where our training and validation data represents.

Convert train and validation AOIs to slippy map tile polygons with supermercado and geopandas

[link to Colab section]

In this step, we’ll use supermercado to generate square tile polygons representing all the slippy map tiles at a specified zoom level that overlap the geoJSON training and validation AOIs we created above.

For this tutorial, we’ll work with slippy map tiles of tile_size=256 and zoom_level=19which yields a manageable number of tiles and satisfactory segmentation results without too much preprocessing or model training time.

You could also try setting a higher or lower zoom_level which would generate more or less tiles at higher or lower resolutions respectively.

Here is an example of different tile zoom_levels over the same area of Zanzibar (see the round, white satellite TV dish for a consistently sized visual reference):

Learn more about slippy maps here, here, and here.

Then we’ll merge our supermercado-generated slippy map tile polygons into a GeoDataFrame with geopandas. We’ll also check for and reconcile overlapping train and validation tiles which would otherwise throw off how we evaluate our progress with model training.

Load slippy map tile image from COG with rio-tiler and corresponding label with geopandas

[link to Colab section]

Now we’ll use rio-tiler and the slippy map tile polygons generated by supermercado to test load a single 256x256 pixel tile from our znz001 COG image file. We will also load the znz001 geoJSON labels into a geopandas GeoDataFrame and crop the building geometries to only those that intersect the bounds of the tile image:

Here is a great introduction to COGs, rio-tiler, and exciting developments in the cloud-native geospatial toolbox by Vincent Sarago of Development Seed:

COG Talk — Part 1: What’s new?

We’ll then create our corresponding 3-channel RGB mask by passing these cropped geometries to solaris’ df_to_px_mask function. Pixel value of 255 in the generated mask:

• in the 1st (Red) channel represent building footprints,
• in the 2nd (Green) channel represent building boundaries (visually looks yellow on the RGB mask display because the pixels overlap red and green+red=yellow),
• and in the 3rd (Blue) channel represent close contact points between adjacent buildings

Make and save all the image and mask tiles

[link to Colab section]

Now that we’ve successfully loaded one tile image from COG with rio-tiler and created its 3-channel RGB mask with solaris, let’s generate our full training and validation datasets. We’ll write some functions and loops to run through all of our trnand val tiles at zoom_level=19 and save them as lossless png files in the appropriate folders with a filename schema of {save_path}/{prefix}{z}_{x}_{y} so we can easily identify and geolocate what tile each file represents.

Train u-net segmentation model with fastai & pytorch

As our deep learning framework and library of tools, we’ll use the excellent fastai library built on top of PyTorch. For more info:

• about Fast.ai, the organization: https://www.fast.ai/about/
• direct links to the free MOOC series:
  - Part 1 (“Practical Deep Learning for Coders”): https://course.fast.ai/index.html
  - Part 2 (“Deep Learning from the Foundations”): https://course.fast.ai/part2

Download and install fastai, set Colab to GPU runtime if needed

[link to Colab section]

Let’s download, install, and set up fastai v1 (currently at 1.0.55). And if we’re not already on it, let’s reset Colab to a GPU runtime (this removes locally stored files since it switches to a new environment so you will have to re-download and untar the training dataset created in above steps):

Colab’s free GPUs range from a Tesla K80, T4, or T8 depending on their availability. See the ===Hardware=== section of show_install() for what GPU type and how much GPU memory is available which will affect the batch size and training time.

For all of these GPUs and mem sizes, a batch size of bs=16 at size=256 should train at <2 mins/epoch without encountering out-of-memory issues but if it does comes up, lower the bs to 8.

Set up data

[link to Colab section]

Now we’ll set up our training dataset of tile images and masks created above to load correctly into fastai for training and validation. The code in this step tracks closely with that of fastai course’s lesson3-camvid so please refer to that lesson video and notebook for more detailed and excellent explanation by Jeremy Howard about the code and fastai’s Data Block API.

The main departures from the camvid lesson notebook is the use of filename string parsing to determine which image and mask files comprise the validation data:

# define the valdation set by fn prefix
holdout_grids = ['znz001val_']
valid_idx = [i for i,o in enumerate(fnames) if any(c in str(o) for c in holdout_grids)]

And we’ll subclass SegmentationLabelList to alter the behavior of open_mask and PIL.Image underlying it in order to open the 3-channel target masks as RGB images (convert_mode=’RGB’) instead of default greyscale 1-channel images (convert_mode=’L’).

We’ll also visually confirm that the image files and channels of the respective target mask file are loaded and paired correctly with a display function show_3ch:

Define custom losses and metrics to handle 3-channel target masks

[link to Colab section]

Here we implement some new loss functions like Dice Loss and Focal Loss which have been shown to perform well in image segmentation tasks. Then we’ll create a MultiChComboLoss class to combine multiple loss functions and calculate them across the 3 channels with adjustable weighting.

The approach of combining a Dice or Jaccard loss to consider image-wide context with individual pixel-focused Binary Cross Entropy or Focal loss with adjustable weighing of the 3 target mask channels has been shown to consistently outperform single loss functions. This is well-documented by Nick Weir’s deep dive into the recent SpaceNet 4 Off-Nadir Building Detection top-5 results:

A deep dive into the SpaceNet 4 winning algorithms

We’ll also adapt our model evaluation metrics (accuracy and dice score) to calculate either a mean score across all channels or for a specified individual channel.

Set up model

[link to Colab section]

We’ll set up fastai’s Dynamic Unet model with an ImageNet-pretrained resnet34 encoder. This architecture, inspired by the original U-net, uses by default many advanced deep learning techniques such as:

• One cycle learning schedule: https://sgugger.github.io/the-1cycle-policy.html
• AdamW optimizer: https://www.fast.ai/2018/07/02/adam-weight-decay/
• Pixel shuffle upsampling with ICNR initiation from super resolution research: https://medium.com/@hirotoschwert/introduction-to-deep-super-resolution-c052d84ce8cf
• Optionally set leaky ReLU, blur, self attention: https://docs.fast.ai/vision.models.unet.html#DynamicUnet

We’ll define our MultiChComboLoss function as a balanced combination of Focal Loss and Dice Loss and set our accuracy and dice metrics:

# set up metrics
acc_ch0 = partial(acc_thresh_multich, one_ch=0)
dice_ch0 = partial(dice_multich, one_ch=0)
metrics = [acc_thresh_multich, dice_multich, acc_ch0, dice_ch0]

# combo Focal + Dice loss with equal channel wts
learn = unet_learner(data, models.resnet34,
                     model_dir='../../models',
                     metrics=metrics, 
                     loss_func=MultiChComboLoss(
                        reduction='mean',
                        loss_funcs=[FocalLoss(gamma=1, alpha=0.95),
                                    DiceLoss()], 
                        loss_wts=[1,1],
                        ch_wts=[1,1,1])
                    )

Also note that our metrics displayed during training shows channel-0 (building footprint channel only) accuracy and dice metrics in the right-most 2 columns while the first two accuracy and dice metrics (left-hand columns) show the mean of the respective metric across all 3 channels.

Train model, inspect results, unfreeze & train more, export for inference

[link to Colab section]

First, we’ll fine-tune our Unet on the decoder part only (leaving the weights for the ImageNet-pretrained resnet34 encoder frozen) for some epochs. Then we’ll unfreeze all the trainable weights/layers of our model and train for some more epochs.

We’ll track the valid_loss, acc_..., and dice_... metrics per epoch as training progresses to make sure they continue to improve and we’re not overfitting. And we set a SaveModelCallback which will track the channel-0 dice score, save a model checkpoint each time there’s an improvement, and reload the highest performing model checkpoint file at the end of training.

We’ll also inspect our model’s results by setting learn.model.eval(), generating some batches of predictions on the validation set, calculating and reshaping the image-wise loss values, and sorting by highest loss first to see the worst performing results (as measured by the loss which may differ in surprising ways from visually gauging results).

Pro-tip: display and view your results every chance you get! You’ll pick up on all kinds of interesting clues about your model’s behavior and how to make it better.

Finally, we’ll export our trained Unet segmentation model for inference purposes as a .pkl file. Learn more about exporting fastai models for inference in this tutorial:

Inference Learner | fastai

Inference on new imagery

[link to Colab section]

With our segmentation model trained and exported for inference use, we will now re-load it as an inference-only model to test on new unseen imagery. We’ll test the generalizability of our trained segmentation model on tiles from drone imagery captured over another part of Zanzibar and in other parts of the world as well as at varying zoom_levels (locations and zoom levels indicated):

We’ll also compare our model inference time per tile on GPU versus CPU:

Post-processing

Predict on a tile, threshold, polygonize, and georegister

[link to Colab section]

For good evaluation of model performance against ground truth, we’ll use another set of labeled data that the model was not trained on. We’ll get this from the larger Zanzibar dataset. Preview the imagery and ground truth labels for znz029 in the STAC browser here:

For demonstration, we’ll use this particular tile at z=19, x=319454, y=270706 from znz029:

Using solaris and geopandas, we’ll convert our model’s prediction as a 3-channel pixel raster output into a GeoJSON file by:

1. thresholding and combining the 3-channels of pixel values in our raw prediction output into a 1 channel binary pixel mask
2. polygonizing this binary pixel mask into shape vectors representing the predicted footprint of every building
3. georegistering the x, y display coordinates of these vectorized building shapes into longitude, latitude coordinates

Evaluate prediction against ground truth

[link to Colab section]

Finally with georegistered building predictions as a GeoJSON file, we can evaluate it against the ground truth GeoJSON file for the same tile.

We’ll clip the ground truth labels to the bounds of this particular tile and use solaris’s Evaluator to calculate the precision, recall, and F1 score. We will also visualize our predicted buildings (in red) against the ground truth buildings (in blue) in this particular tile:

For more information about these common evaluation metrics for models applied to overhead imagery, see the following articles and more by the SpaceNet team:

• The SpaceNet Metric
• The good and the bad in the SpaceNet Off-Nadir Building Footprint Extraction Challenge

Ideas to Try for Performance Gains

Congratulations, you did it!

You’ve completed the tutorial and now know how to do everything from producing training data to creating a deep learning model for segmentation to postprocessing and evaluating your model’s performance.

To flex your newfound knowledge and make your model perform potentially much better, try implementing some or all these ideas:

• Create and use more training data: there are 13 grids’ worth of training data for Zanzibar released as part of the Open AI Tanzania Building Footprint Segmentation Challenge dataset.
• Change the zoom_level of your training/validation tiles. Better yet, try using tiles across multiple zooms (i.e. z21, z20, z19, z18). Note that with multiple zoom levels over the same imagery, you should be extra careful of overlapping tiles across those different zoom levels. ← test your understanding of slippy map tiles by checking that you understand what I mean here but feel free to message me for the answer!
• Change the Unet’s encoder to a bigger or different architecture (i.e. resnet50, resnet101, densenet).
• Change the combinations, weighting, and hyperparameters of the loss functions. Or implement completely new loss functions like Lovasz Loss.
• Try different data augmentation combinations and techniques.
• Train for more epochs and with different learning rate schedules. Try mixed-precision for faster model training.
• Your idea here.

I look forward to seeing what you discover!

Coming Up

If you liked this tutorial, look forward to next ones which will potentially cover topics like:

• classifying building completeness (foundation, incomplete, complete)
• inference on multiple tiles and much larger images
• working with messy, sparse, imperfect training data
• model deployment and inference at scale
• examining data/model biases, considerations of fairness, accountability, transparency, and ethics

Curious about more geospatial deep learning topics? Did I miss something? Share your questions and thoughts in the comments so I can add them into this and next tutorials.

Good luck and happy deep learning!

Acknowledgments and Special Thanks to

• World Bank GFDRR’s Open Data for Resilience Initiative (OpenDRI) for consultation projects which have inspired & informed.
• Zanzibar Mapping Initiative, OpenAerialMap, State University of Zanzibar (SUZA), Govt of Zanzibar’s Commission for Lands, & WeRobotics for the 2018 Open AI Tanzania Building Footprint Segmentation Challenge.
• Fast.ai’s team, contributors, & community for both “making neural nets uncool again” and pushing its cutting edge (very cool).
• SpaceNet & Cosmiq Works for the open challenges, datasets, knowledge-sharing, Solaris geoML toolkit, & more that advance geospatial machine learning.
• Contributors to COG, STAC, and more initiatives advancing the cloud native geospatial ecosystem.
• Free & open source creators & collaborators everywhere for the invaluable public goods you provide.

Imagine

You open Google Maps and enter “coffee” to find shops nearby. The app proceeds to download a map of your entire city at the highest detail. You wait minutes and 100s of MBs download to your phone before 4 or 5 location pins drop closest to you.

If the app did this every time you’re in a new area or search for something different, you would probably stop using it. This scenario is extreme, even absurd, yet we often do something similar with geospatial data in deep learning.

We download full-sized satellite or aerial imagery (at 100s of MBs to GBs per image or per band), crop, resize, & tile them to the areas, sizes, & formats we need, and run our model training or inference on the end product while holding a relatively large portion of the source data unused.

If we need the highest resolution & full coverage of all images to train our models, this works. But what if we want to evaluate models on small select subareas of new images or analyze very large areas at faster speed & less detail than what the original file provides?

Could we access just the relevant areas at only the resolutions we need? Until recently, we had little choice but to download the entirety of every image. Now, thanks to Cloud Optimized GeoTIFFs (COG), we have a better way to run deep learning models on geospatial data more efficiently at any size & scale.

A Brief Intro to COG

In the geospatial world, many exciting developments are moving us towards more interoperable, cloud-native architectures for data processing & analysis:

• Cloud-Optimized GeoTIFF (COG)
• Spatial Temporal Asset Catalog (STAC)
• Analysis Ready Data (ARD)

The whys, whats, and hows of these related initiatives have been well explained by Chris Holmes in his 3-part “Cloud-Native Geoprocessing” series so I won’t get into that here. I will touch on these as they pertain to geospatial deep learning in subsequent posts.

Today, though, is all about the COG. Summed up in an excellent presentation on efficient cloud workflows by @echeipesh from GeoTrellis, we want this:

“Hey, let’s not download the whole file every time.”

COGs deliver: accessing huge geospatial files becomes a speedy & selective data-streaming experience using the same web tech that enable videos to start playing before the whole file is downloaded. Here are some recent advances & implementations showing what’s newly possible:

COG for Deep Learning

The functionality we’ll focus on for deep learning is how COGs use tiling & overviews:

Tiling creates a number of internal ‘tiles’ inside the actual image, instead of using simple ‘stripes’ of data. With a stripe of data then the whole file needs to be read to get the key piece. With tiles much quicker access to a certain area is possible, so that just the portion of the file that needs to be read is accessed.

Overviews create downsampled versions of the same image. This means it’s ‘zoomed out’ from the original image — it has much less detail (1 pixel where the original might have 100 or 1000 pixels), but is also much smaller. Often a single GeoTIFF will have many overviews, to match different zoom levels. These add size to the overall file, but are able to be served much faster, since the renderer just has to return the values in the overview instead of figuring out how to represent 1000 different pixels as one.

The organization of these tiles & overviews delivered by a COG tile server (like tiles.rdnt.io thanks to Radiant.Earth & Seth Fitzsimmons) generally follows the slippy map tile naming convention:

• Tiles are 256 × 256 px or 512 × 512 px PNG or JPG files
• Filename(url) format is /{zoom}/{x}/{y}.png
• Each zoom level is a directory {zoom}, each column is a subdirectory {x}, and each tile in that column is a file {y}

Hmm, 512 px & 256 px square images. Starting to sound familiar?

COG tiles & overviews hand us geospatial data on a platter: consistently georeferenced & internally organized, optimized for fast access & visibility at every zoom level, and formatted in a familiar way for deep learning models.

To test drive the potential, I created an Input COG → Model Inference → Output COG workflow that:

1. gets overview tiles from any COG at any zoom level (or multiple levels)
2. runs inference on each tile and reassembles the results
3. saves the output as a properly georeferenced and validated COG file

An Example, Step-by-Step

Here’s that workflow in action. We’ll use this 7-cm resolution drone image taken by the UN’s International Org. for Migration of the Rohingya refugee camps near Cox’s Bazar, Bangladesh hosted as a COG on OpenAerialMap:

1. Get tiles at zoom levels 17, 18, & 19 served by tiles.rdnt.io

2. Run inference per tile with model trained offline to find built-up areas* (binary semantic segmentation)

* model used for demo of workflow; not fitted to this data so results may appear suboptimal. Model training & inference will be covered in a later post.

3. Reassemble tiled results into full output map at each zoom level

4. Ensemble into final output map (with new color range)

5. Calculate geo bounds & save as georeferenced COG

How Fast Is It?

Running the example workflow on a remote GPU instance (Paperspace’s P5000 machine) took ~30 seconds:

starting inference for zoom level: 17

100%|██████████| 15/15 [00:01<00:00,  9.36it/s]

starting inference for zoom level: 18

100%|██████████| 50/50 [00:10<00:00,  4.78it/s]

starting inference for zoom level: 19

100%|██████████| 152/152 [00:13<00:00, 11.04it/s]

CPU times: user 14.2 s, sys: 968 ms, total: 15.1 s
Wall time: 25.9 s

I was also impressed at the reasonable speed using CPU only: ~12 minutes on my run-of-mill 2015 Macbook Pro even though code is not optimized for performance:

CPU times: user 14min 51s, sys: 1min 26s, total: 16min 17s
Wall time: 11min 43s

It’s faster because using COG overview tiles at specific resolutions (zoom levels) directly gets the information we need while avoiding the unnecessary & heavy data management typical of preparing geospatial data for deep learning.

In the 1st example, the original COG file size is only 30MB which wouldn’t be a very perceptible difference. Using a much larger example, this 900MB source image of a bigger area was processed in 90 seconds on GPU. Working traditionally, I would still be downloading the file (at typical broadband speeds):

A zoomed-in & overlaid view of the same output as above shows actionable details (where built-up areas are) are preserved up close:

Oh, the Possibilities

COGs enable us to run our deep learning models more rapidly, lightly, & simply on geospatial data at any size or scale.

This data includes satellite & aerial imagery from Landsat on AWS, OpenAerialMap, & Planet and more soon as providers increasingly go to COG:

The advance of COG for deep learning means that we could:

• work selectively with any-sized subareas and as many zoom levels of source imagery as we need to get useful results.
• select single spectral bands or mix-and-match any band combination with one change in the tile server parameter (i.e. rgb=1,1,1).
• create new models by serving COG tiles directly into our training pipelines with labels generated on the fly, perhaps via geospatial machine learning data prep tools such as Robosat or Label Maker.
• test many models on imagery of one area, or one model on many areas and a wide gamut of visual conditions to evaluate their generalizability to real world data.
• deploy models for any COG data provider with less cost & infrastructure (cloud-based, CPU only), making it more feasible to AI-enhance many localized humanitarian, environmental, & community-based geospatial projects like those being carried out by WeRobotics’ Flying Labs.

In following posts, we’ll cover this workflow in technical detail (with code examples), experiment with these possibilities, and encourage more new ideas for cloud-native geospatial deep learning.

I look forward to seeing & sharing what you come up with!

Like what you’re reading? Want to protect our health & prepare our communities for climate change? If you’re looking to do your best work in geospatial analysis & deep learning to tackle our hardest systems challenges in environmental health & justice, Anthropocene Labs is looking for you! Get in touch with dave(at)anthropo dot co

The fire burned for more than a month, though its spread was contained several weeks ago. Heavy rains earlier this week, which caused land burned by the fire to create mudflows that buried neighborhoods, helped fully extinguish the blaze. In the end, the fire burned 281,893 acres.

The fire eclipsed the 2003 Cedar fire in San Diego County, which burned 273,246 acres.

The milestone reaffirmed 2017 as the most destructive fire season in the state. In October, a series of fires in wine country burned more than 10,000 homes and killed more than 40 people.

Those blazes, along with the Thomas fire, were fueled by dry conditions and intense winds.

http://www.latimes.com/local/lanow/la-me-thomas-fire-contained-20180112-story.html

License, Sources, and Technical Notes:

• Animated GIF and Juxtapose images licensed under a Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License
• Landsat-8 bands combination: 7 (shortwave infrared 2), 5 (near infrared), 2 (blue)
• Source GeoTiff files from Google Cloud Public Datasets:
  ‘LC08_L1TP_042036_20171022_20171107_01_T1’,
  ‘LC08_L1TP_042036_20171209_20171223_01_T1’,
  ‘LC08_L1TP_042036_20171225_20180103_01_T1’,
  ‘LC08_L1TP_042036_20180126_20180207_01_T1’
• Made in Python with Jupyter notebooks, sat-utils, geojson.io, rasterio, ezgif.com, juxtaposeJS

281,893 Acres was originally published in Delta Anthropoco on Medium, where people are continuing the conversation by highlighting and responding to this story.

Levees stopped the natural flow of fresh water and sediment that reinforced the fragile marshes. Oil and gas companies dredged through the mud to lay pipelines and build canals, carving paths for saltwater to intrude and kill the freshwater vegetation that held the land together. The unstoppable, glacial momentum of sea-level rise has only made things worse. Today, almost nothing remains of what was very recently a vast expanse of bountiful marshes and swampland.

Isle de Jean Charles, home to the Biloxi-Chitimacha-Choctaw band of Native Americans, has lost 98 percent of its land since 1955. Its 99 remaining residents have been dubbed “America’s first climate refugees.”

…

The residents of Isle de Jean Charles won’t be alone in their exodus. There will be up to 13 million climate refugees in the United States by the end of this century. Even if humanity were to stop all carbon emissions today, at least 414 towns, villages, and cities across the country would face relocation, according to a study published in the Proceedings of the National Academy of Sciences. If the West Antarctic Ice Sheet collapses, researchers predict that the number will exceed 1,000.

And this isn’t a distant threat. At least 17 communities, most of which are Native American or Native Alaskan, are already in the process of climate-related relocations. Yet despite its inevitability, there is no official framework to handle this displacement. There is no U.S. government agency, process, or funding dedicated to confronting this impending humanitarian crisis.

Only one climate-related relocation is currently funded and administered by the government: the Isle de Jean Charles Resettlement Project.

- “How to Save a Town From Rising Waters” [Michael Isaac Stein at CityLab]

Learn more at:

• https://vimeo.com/ondemand/cantstopthewater
• http://www.coastalresettlement.org/
• https://www.motherjones.com/environment/2017/10/climate-refugees-trump-hud/

License, Sources, and Technical Notes:

• Animated GIF at top licensed under a Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License
• Landsat-8 bands combination: 5 (near infrared), 2 (blue), 1 (coastal blue)
• Source GeoTiff files from Google Cloud Public Datasets:
  ‘LC08_L1TP_022040_20131218_20170307_01_T1’,
  ‘LC08_L1TP_022040_20140409_20170307_01_T1’,
  ‘LC08_L1TP_022040_20141119_20170302_01_T1’,
  ‘LC08_L1TP_022040_20150207_20170301_01_T1’,
  ‘LC08_L1TP_022040_20150919_20170225_01_T1’,
  ‘LC08_L1TP_022040_20160210_20170224_01_T1’,
  ‘LC08_L1TP_022040_20160226_20170224_01_T1’,
  ‘LC08_L1TP_022040_20160313_20170224_01_T1’,
  ‘LC08_L1TP_022040_20161210_20170219_01_T1’,
  ‘LC08_L1TP_022040_20171026_20171107_01_T1’,
  ‘LC08_L1TP_022040_20171127_20171206_01_T1’,
  ‘LC08_L1TP_022040_20180114_20180120_01_T1’,
• Made in Python with Jupyter notebooks, sat-utils, geojson.io, rasterio, ezgif.com

Then we had water, all the time. was originally published in Delta Anthropoco on Medium, where people are continuing the conversation by highlighting and responding to this story.

Time-lapse satellite imagery of two agricultural Californian counties over 2017 in natural color.

Left: Imperial County in southeastern California next to the Salton Sea.

Right: Sutter County in northern California’s Sacramento Valley.

A contrast of landscape, geography, and childhood asthma rates:

License, Sources, and Technical Notes:

• First two animated GIFs licensed under a Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License
• Landsat-8 bands combination: 4 (red), 3 (green), 2 (blue)
• Source GeoTiff files from Google Cloud Public Datasets:
  Imperial:
  ‘LC08_L1TP_039037_20170118_20170218_01_T1’,
  ‘LC08_L1TP_039037_20170307_20170317_01_T1’,
  ‘LC08_L1TP_039037_20170526_20170615_01_T1’,
  ‘LC08_L1TP_039037_20170729_20170811_01_T1’,
  ‘LC08_L1TP_039037_20170915_20170928_01_T1’,
  ‘LC08_L1TP_039037_20171118_20171205_01_T1’,
  Sutter:
  ‘LC08_L1TP_044033_20170121_20170218_01_T1’,
  ‘LC08_L1TP_044033_20170427_20170515_01_T1’,
  ‘LC08_L1TP_044033_20170614_20170628_01_T1’,
  ‘LC08_L1TP_044033_20170817_20170826_01_T1’,
  ‘LC08_L1TP_044033_20171004_20171014_01_T1’,
  ‘LC08_L1TP_044033_20171207_20171223_01_T1’,
• Made in Python with Jupyter notebooks, geojson.io, rasterio, ezgif.com

California Breathing, 2017 in Two Counties was originally published in Delta Anthropoco on Medium, where people are continuing the conversation by highlighting and responding to this story.

Annotated version:

A time-lapse of satellite images showing Cape Town’s largest fresh water reservoir (Theewaterskloof Dam) depleting to historical lows after 3 successive dry years.

Images are generated as a false color combination of Near Infrared, Green, and Coastal/Aerosol (Landsat-8 bands 5, 3, 1) to highlight vegetation (red, orange) and coastal/inland waters (shades of blue).

As of Feb 2, 2018, the reservoir is at 12.2% capacity. The city of Cape Town currently forecasts Day Zero — “the day the taps will be turned off” — to arrive Apr 16 2018 . If and when that day comes, the city plans to provide 200 water distribution stations to allocate 25 liters per person per day. The average American uses over 300 liters of water per day.

Learn more below:

“Cape Town’s Reservoirs Are Getting Terrifyingly Low”

[Earther.com]

“As things stand, the challenge exceeds anything a major city has had to face anywhere in the world since the Second World War or 9/11,” Helen Zille, the Premier of the Western Cape, wrote in an op-ed late last month. “I personally doubt whether it is possible for a city the size of Cape Town to distribute sufficient water to its residents, using its own resources, once the underground waterpipe network has been shut down.”

The city got to this point because rains have failed the past three years. Last year was Cape Town’s driest year on record, taking the record from 2016. Oh, and 2015 was the fourth-driest year ever recorded at the city’s main weather station at the airport.

“How severe is this drought, really?”

[Piotr Wolski at Climate System Analysis Group]

Some time ago, I’ve heard a story told by an old Kenyan herdsman. The Madala said that every now and then, his village experienced a cattle-killing drought. Every decade, or so, they experienced a goat-killing drought (goats being more sturdy than cows, would die only if conditions were more extreme). And once in a man’s (woman’s of course too) lifetime, they would experience a man-killing drought.

So is Cape’s drought cattle-, goat- or man-killing one? This question can, of course, be rephrased in scientific terms — what is the return interval of the drought we are experiencing? Which basically means: how often, on average, can we expect a drought of a magnitude of the one that we have now, or more severe, to occur?

“Cape Town’s Water is Running Out“

[NASA Earth Observatory Image of the Day]

Piotr Wolski, a hydrologist at the Climate Systems Analysis Group at the University of Cape Town, has analyzed rainfall records dating back to 1923 to get a sense of the severity of the current drought compared to historical norms. His conclusion is that back-to-back years of such weak rainfall (like 2016–17) typically happens about once just every 1,000 years.

Population growth and a lack of new infrastructure has exacerbated the current water shortage. Between 1995 and 2018, the Cape Town’s population swelled by roughly 80 percent. During the same period, dam storage increased by just 15 percent.

The city did recently accelerate development of a plan to increase capacity at Voëlvlei Dam by diverting winter rainfall from the Berg River. The project had been scheduled for completion in 2024, but planners are now targeting 2019. The city is also working to build a series of desalination plants and to drill new groundwater wells that could produce additional water.

License, Sources, and Technical Notes:

• First two animated GIFs licensed under a Creative Commons Attribution-NonCommercial-ShareAlike 4.0 International License
• Landsat-8 bands combination: 5 (near infrared), 3 (green), 1 (coastal/aerosol)
• Source GeoTiff files from Google Cloud Public Datasets:
  ‘LC08_L1TP_175084_20131218_20170427_01_T1’,
  ‘LC08_L1TP_175084_20141018_20170418_01_T1’,
  ‘LC08_L1TP_175084_20150223_20170412_01_T1’,
  ‘LC08_L1TP_175084_20160109_20170405_01_T1’,
  ‘LC08_L1TP_175084_20161226_20170315_01_T1’,
  ‘LC08_L1TP_175084_20171010_20171024_01_T1’,
  ‘LC08_L1TP_175084_20171127_20171206_01_T1’,
  ‘LC08_L1TP_175084_20180114_20180120_01_T1’
• Made in Python with Jupyter notebooks, geojson.io, rasterio, ezgif.com

Once Every 1,000 Years was originally published in Delta Anthropoco on Medium, where people are continuing the conversation by highlighting and responding to this story.

“Two billion years ago, cyanobacteria oxygenated the atmosphere and powerfully disrupted life on Earth. But they didn’t know it.

We’re the first species that’s become a planet-scale influence and is aware of that reality. That’s what distinguishes us.”

- Andrew Revkin (@revkin) via Smithsonian

Exploring anthropos’ influence and raising our fledging self-awareness, starting with this one human. Hi[at]anthropo[dot]co.