Fluvial erosion linked to warming in the Arctic

This repository contains code related to Temperature Controls Erosion along Arctic Rivers by Jordan Fields, Jack Kreisler, Evan Dethier, John Perrotti, and Carl Renshaw.

Datasets

The primary dataset for this study was Landsat 5, 7, and 8 surface reflectance images (Collection 2). Additional datasets include:

- MERIT Hydro
- Arctic DEM
- SRTM
- HydroLAKES
- ERA5-Land daily
- Downscaled CMIP6 GCM 

Software

Datasets were accessed and analyzed through Google Earth Engine. Results were exported from Earth Engine and further analyzed in R for trend analysis, statistics, and visualizations. We used the data.table package for it's speed during our data organization but created plots using ggplot. 

Code

All code needed to reproduce the results is provided here. The purpose of and any results generated by each file are briefly described here.

1. Google Earth Engine

buffer_generator

Generates centroids along the Mackenzie River stream network inside the project study area. Exports MERIT_Segment_Buffers (which contains the centroids) and MERIT_Unioned_Streams (used in thiessen_polys to buffer the entire stream network). The buffers themselves are not used at all, only the centroids.

thiessen_polys

Generates Thiessen polygons along the entire stream network within the study area. These polygons allow us to define discrete, non-overlapping study sites of generally similar size without any gaps between sites. Exports the collection of Thiessen polygons.

That collection has been manually reviewed to remove any problematic polygons, such as those that cover lakes, those with no water area, and those in areas where the stream is to narrow to consistently detect water (~100 m). The edited collection, used for all analysis, is Good_Thiessens_3.

cpixel_generator

Generates a collection of yearly images that allow us to mask out and ignore water in near-channel lakes within New_cpixel_AOI. This region covers the study area north of ~67.5 degrees of latitude, where near-channel lakes are particularly numerous. Each year's image is created with a connected pixels analysis where any water pixels with fewer than 1024 connected neighbors are masked out and therefore not taken into account for any analysis. Exports cpixel_Masks, the collection of yearly masks.

Note -- review of yearly masks showed significant variation in the connectivity of near-channel lakes, which led to a noisy land gained/lost signal in the cpixel AOI. Thus, results from within the cpixel AOI are only analyzed in bulk analysis of all results and not included in final reported findings in the main text (as noted in Online Methods).

slope_mask

Creates a single image to mask out all areas with a slope greater than 22 within High_Relief_Region (an asset available in the script). Manual review of water detection accuracy showed that mountain shadows were sometimes erroneously detected as water, so creating and applying this mask allows for more accurate water detection. Exports Mack_SlopeMask.

Groundtruthing

Used to manually delineate water for selected ground truth sites, with results shown Online Methods Fig. 3

MASTER_Temp_Calcs

Uses ERA5 reanalysis data to calculate degree days and days above 0°C for each year for each Thiessen polygon. Exports the mean value of degree days and day above 0°C within each polygon in each 5 year window (i.e. every thiessen polygon end up with values for every year from 1990-2023). Exports these values to Google Drive in batches, as running all at once exceeds Earth Engine memory limits. The files are aggregated in R (data.table) and the data included in the file: Master1_ArcticErosion.csv

MASTER_FUTURE_Temp_Calcs

Uses CMIP6 model temperature data and the same methods as "Master_Temp_Calcs" to produce mean value of degree days and day above 0°C within each polygon in each 5 year window (i.e. every thiessen polygon end up with values for every year from 1990-2023). Exports these values to Google Drive in batches, as running all at once exceeds Earth Engine memory limits. The files are aggregated in R (data.table) and compiled in the file: Future_Temperature.csv

MASTER_NDVI_Calcs

Calculates mean per-pixel NDVI values for each year, then computes the mean NDVI value within each Thiessen polygon for each year. Exports these values to Google Drive in batches, as running all at once exceeds Earth Engine memory limits. The files are aggregated in R (data.table) and compiled into the file: NDVI_all_1.csv

MASTER_meander_rates

The main analysis for this study. Detects water in all Landsat 5, 7, and 8 images within the study area using a Water Ratio Index (WRI) algorithm supplemented by band-based masks, the connected pixel masks, and the slope mask. Determines yearly water occurrence across the entire study area using a two-year image window, where water pixels are those that are water in at least 20% of the images. Calculates area that changed from land to water (or vice versa) over rolling 5-year windows from 1990-2023. Counts total water area, land gained, and land lost within each Thiessen polygon for each rolling 5-year window. Exports these values to Google Drive in batches, as running all at once exceeds Earth Engine memory limits. The files are aggregated in R (data.table) and and the data (raw values and discharge-corrected values) included in the file: Master1_ArcticErosion.csv 

Assets

To reduce computational complexity and prevent time-out errors, some items (such as Thiessen polygons and connected pixel masks) are computed ahead of time and imported into Earth Engine scripts. These assets are listed in list_of_assets, are publicly available on Google Earth Engine as Cloud Assets, thereby accessible by users via the GEE scripts included here, but are also inlcuded in the folder exported_assets.zip above so that users can upload them directly, if neeed. 

Batch Exports

Due to the vast number of images used and despite the ability to run computations on Earth Engine servers, we still encountered some memory limitations. To avoid these limitations and allow results to compute, in some cases we export the NDVI/temperature/meander results in batches. Each batch export comprises all years but only a subset of all Thiessen polygons. The for-loops and other commands used for batch exports are present in scripts (for transparency) but commented out (for clarity). Reproducing our results may or may not require batch exports, depending on code optimizations.

2. R Scripts

Mack_River_DataAnalysis_CLEAN.R
This script contains all the data organization and analysis done outside of Google Earth Engine (GEE). A few of the critical actions that are done with these scripts are:

• Aggregation of all the batch exported files from GEE
• Discharge corrections to area lost values
• 5yr rolling values for NDVI values
• Summarization of data by year, latitude, stream order and permafrost extent
• Plotting and statistical analyses

The inputs used in this script are included in the folder: Input_R_Data.zip. These inputs include:

• Raw GEE Exports files in folder GEE_Exports
• Discharge Data in folder Discharge_Data
• Suspended Sediment Data mackenzie_river_ssc_warm_months_summary.csv
• Permafrost Data Thiessen_Permafrost_Data_0926.csv
• Groundtruth Data Groundtruth_Data.csv
• Case study data Manually_Measured_Data.csv

The analysis in R creates the following outputs, which are included in the folder: R_Exports.zip

• All data Master1_ArcticErosion.csv
• Future temperature data fro CMIP6 Future_Temperature.csv
• Data summarized by year Summary1_ArcticErosion_Year.csv
• Data summarized by stream order Summary2_ArcticErosion_StrmOrd.csv
• Data summarized by latitude Summary3_ArcticErosion_LatBins.csv
• Data summarized by permafrost extent Summary4_ArcticErosion_Permafrost_Extent.csv
• Data summarized by thaw depth Summary5_ArcticErosion_Thaw_Depth.csv
• Total water area in each thiessen Total_Area_Oct25_23.csv