diff --git a/scripts/run.sh b/scripts/run.sh index cfc8f75..3d6c9b7 100755 --- a/scripts/run.sh +++ b/scripts/run.sh @@ -22,6 +22,7 @@ fi declare -a SCENARIOS=("arable" "restore" "restore_all" "urban" "pasture" "restore_agriculture") declare -a TAXAS=("AMPHIBIA" "AVES" "MAMMALIA" "REPTILIA") +export CURVE=0.25 python3 ./prepare_layers/generate_crosswalk.py --output "${DATADIR}"/crosswalk.csv @@ -162,7 +163,7 @@ done # Generate the batch job input CSVs python3 ./utils/speciesgenerator.py --datadir "${DATADIR}" --output "${DATADIR}"/aohbatch.csv -python3 ./utils/persistencegenerator.py --datadir "${DATADIR}" --output "${DATADIR}"/persistencebatch.csv +python3 ./utils/persistencegenerator.py --datadir "${DATADIR}" --curve "${CURVE}" --output "${DATADIR}"/persistencebatch.csv # Calculate all the AoHs littlejohn -j 700 -o "${DATADIR}"/aohbatch.log -c "${DATADIR}"/aohbatch.csv "${VIRTUAL_ENV}"/bin/python3 -- ./aoh-calculator/aohcalc.py --force-habitat @@ -189,14 +190,14 @@ for SCENARIO in "${SCENARIOS[@]}" do for TAXA in "${TAXAS[@]}" do - python3 ./utils/raster_sum.py --rasters_directory "${DATADIR}"/deltap/"${SCENARIO}"/0.25/"${TAXA}"/ --output "${DATADIR}"/deltap_sum/"${SCENARIO}"/0.25/"${TAXA}".tif + python3 ./utils/raster_sum.py --rasters_directory "${DATADIR}"/deltap/"${SCENARIO}"/"${CURVE}"/"${TAXA}"/ --output "${DATADIR}"/deltap_sum/"${SCENARIO}"/"${CURVE}"/"${TAXA}".tif done - python3 ./utils/species_totals.py --aohs "${DATADIR}"/deltap/"${SCENARIO}"/0.25/ --output "${DATADIR}"/deltap/"${SCENARIO}"/0.25/totals.csv + python3 ./utils/species_totals.py --aohs "${DATADIR}"/deltap/"${SCENARIO}"/"${CURVE}"/ --output "${DATADIR}"/deltap/"${SCENARIO}"/"${CURVE}"/totals.csv # Generate final map - python3 ./deltap/delta_p_scaled.py --input "${DATADIR}"/deltap_sum/"${SCENARIO}"/0.25/ \ + python3 ./deltap/delta_p_scaled.py --input "${DATADIR}"/deltap_sum/"${SCENARIO}"/"${CURVE}"/ \ --diffmap "${DATADIR}"/habitat/"${SCENARIO}"_diff_area.tif \ - --totals "${DATADIR}"/deltap/"${SCENARIO}"/0.25/totals.csv \ - --output "${DATADIR}"/deltap_final/scaled_"${SCENARIO}"_0.25.tif + --totals "${DATADIR}"/deltap/"${SCENARIO}"/"${CURVE}"/totals.csv \ + --output "${DATADIR}"/deltap_final/scaled_"${SCENARIO}"_"${CURVE}".tif done diff --git a/scripts/slurm.sh b/scripts/slurm.sh new file mode 100644 index 0000000..25d167d --- /dev/null +++ b/scripts/slurm.sh @@ -0,0 +1,207 @@ +#!/bin/bash +# +# Assumes you've set up a python virtual environement in the current directory. +# +# In addition to the Python environemnt, you will need the following extra command line tools: +# +# https://github.com/quantifyearth/reclaimer - used to download inputs from Zenodo directly +# https://github.com/quantifyearth/littlejohn - used to run batch jobs in parallel + +set -e + +source ${HOME}/venvs/life/bin/activate +cd ${HOME}/dev/life +export PATH=$PATH:$HOME/go/bin + +if [ -z "${DATADIR}" ]; then + echo "Please specify $DATADIR" + exit 1 +fi + +if [ -z "${VIRTUAL_ENV}" ]; then + echo "Please specify run in a virtualenv" + exit 1 +fi + +declare -a SCENARIOS=("arable" "restore" "restore_all" "urban" "pasture" "restore_agriculture") +declare -a TAXAS=("AMPHIBIA" "AVES" "MAMMALIA" "REPTILIA") +export CURVE=gompertz + + +python3 ./prepare_layers/generate_crosswalk.py --output "${DATADIR}"/crosswalk.csv + +# Get habitat layer and prepare for use +reclaimer zenodo --zenodo_id 4058819 \ + --filename iucn_habitatclassification_composite_lvl2_ver004.zip \ + --extract \ + --output "${DATADIR}"/habitat/jung_l2_raw.tif + +reclaimer zenodo --zenodo_id 4058819 \ + --filename lvl2_changemasks_ver004.zip \ + --extract \ + --output "${DATADIR}"/habitat/ + +python3 ./prepare_layers/make_current_map.py --jung "${DATADIR}"/habitat/jung_l2_raw.tif \ + --update_masks "${DATADIR}"/habitat/lvl2_changemasks_ver004 \ + --crosswalk "${DATADIR}"/crosswalk.csv \ + --output "${DATADIR}"/habitat/current_raw.tif \ + -j 16 + +python3 ./aoh-calculator/habitat_process.py --habitat "${DATADIR}"/habitat/current_raw.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat_maps/current/ + +# Get PNV layer and prepare for use +reclaimer zenodo --zenodo_id 4038749 \ + --filename pnv_lvl1_004.zip \ + --extract \ + --output "${DATADIR}"/habitat/pnv_raw.tif + +python3 ./aoh-calculator/habitat_process.py --habitat "${DATADIR}"/habitat/pnv_raw.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat_maps/pnv/ + + +# Generate an area scaling map +python3 ./prepare_layers/make_area_map.py --scale 0.016666666666667 --output "${DATADIR}"/area-per-pixel.tif + +# Generate the arable scenario map +python3 ./prepare_layers/make_arable_map.py --current "${DATADIR}"/habitat/current_raw.tif \ + --output "${DATADIR}"/habitat/arable.tif + +python3 ./aoh-calculator/habitat_process.py --habitat "${DATADIR}"/habitat/arable.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat_maps/arable/ + +python3 ./prepare_layers/make_diff_map.py --current "${DATADIR}"/habitat/current_raw.tif \ + --scenario "${DATADIR}"/habitat/arable.tif \ + --area "${DATADIR}"/area-per-pixel.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat/arable_diff_area.tif + +# Generate the pasture scenario map +python3 ./prepare_layers/make_pasture_map.py --current "${DATADIR}"/habitat/current_raw.tif \ + --output "${DATADIR}"/habitat/pasture.tif + +python3 ./aoh-calculator/habitat_process.py --habitat "${DATADIR}"/habitat/pasture.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat_maps/pasture/ + +python3 ./prepare_layers/make_diff_map.py --current "${DATADIR}"/habitat/current_raw.tif \ + --scenario "${DATADIR}"/habitat/pasture.tif \ + --area "${DATADIR}"/area-per-pixel.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat/pasture_diff_area.tif + +# Generate the restore map +python3 ./prepare_layers/make_restore_map.py --pnv "${DATADIR}"/habitat/pnv_raw.tif \ + --current "${DATADIR}"/habitat/current_raw.tif \ + --crosswalk "${DATADIR}"/crosswalk.csv \ + --output "${DATADIR}"/habitat/restore.tif + +python3 ./aoh-calculator/habitat_process.py --habitat "${DATADIR}"/habitat/restore.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat_maps/restore/ + +python3 ./prepare_layers/make_diff_map.py --current "${DATADIR}"/habitat/current_raw.tif \ + --scenario "${DATADIR}"/habitat/restore.tif \ + --area "${DATADIR}"/area-per-pixel.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat/restore_diff_area.tif + +# Generate the restore map +python3 ./prepare_layers/make_restore_agriculture_map.py --pnv "${DATADIR}"/habitat/pnv_raw.tif \ + --current "${DATADIR}"/habitat/current_raw.tif \ + --crosswalk "${DATADIR}"/crosswalk.csv \ + --output "${DATADIR}"/habitat/restore_agriculture.tif + +python3 ./aoh-calculator/habitat_process.py --habitat "${DATADIR}"/habitat/restore_agriculture.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat_maps/restore_agriculture/ + +python3 ./prepare_layers/make_diff_map.py --current "${DATADIR}"/habitat/current_raw.tif \ + --scenario "${DATADIR}"/habitat/restore_agriculture.tif \ + --area "${DATADIR}"/area-per-pixel.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat/restore_agriculture_diff_area.tif + +# Generate the restore all map +python3 ./prepare_layers/make_restore_all_map.py --pnv "${DATADIR}"/habitat/pnv_raw.tif \ + --current "${DATADIR}"/habitat/current_raw.tif \ + --crosswalk "${DATADIR}"/crosswalk.csv \ + --output "${DATADIR}"/habitat/restore_all.tif + +python3 ./aoh-calculator/habitat_process.py --habitat "${DATADIR}"/habitat/restore_all.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat_maps/restore_all/ + +python3 ./prepare_layers/make_diff_map.py --current "${DATADIR}"/habitat/current_raw.tif \ + --scenario "${DATADIR}"/habitat/restore_all.tif \ + --area "${DATADIR}"/area-per-pixel.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat/restore_all_diff_area.tif + +# Generate urban all map +python3 ./prepare_layers/make_constant_habitat.py --examplar "${DATADIR}"/habitat_maps/arable/lcc_1401.tif \ + --habitat_code 14.5 \ + --crosswalk "${DATADIR}"/crosswalk.csv \ + --output "${DATADIR}"/habitat_maps/urban + +python3 ./prepare_layers/make_constant_diff_map.py --current "${DATADIR}"/habitat/current_raw.tif \ + --habitat_code 14.5 \ + --crosswalk "${DATADIR}"/crosswalk.csv \ + --area "${DATADIR}"/area-per-pixel.tif \ + --scale 0.016666666666667 \ + --output "${DATADIR}"/habitat/urban_diff_area.tif + +# Fetch and prepare the elevation layers +reclaimer zenodo --zenodo_id 5719984 --filename dem-100m-esri54017.tif --output "${DATADIR}"/elevation.tif +gdalwarp -t_srs EPSG:4326 -tr 0.016666666666667 -0.016666666666667 -r max -co COMPRESS=LZW -wo NUM_THREADS=40 "${DATADIR}"/elevation.tif "${DATADIR}"/elevation-max.tif +gdalwarp -t_srs EPSG:4326 -tr 0.016666666666667 -0.016666666666667 -r min -co COMPRESS=LZW -wo NUM_THREADS=40 "${DATADIR}"/elevation.tif "${DATADIR}"/elevation-min.tif + +# Get species data per taxa from IUCN data +for TAXA in "${TAXAS[@]}" +do + python3 ./prepare_species/extract_species_psql.py --class "${TAXA}" --output "${DATADIR}"/species-info/"${TAXA}"/ --projection "EPSG:4326" +done + +# Generate the batch job input CSVs +python3 ./utils/speciesgenerator.py --datadir "${DATADIR}" --output "${DATADIR}"/aohbatch.csv +python3 ./utils/persistencegenerator.py --datadir "${DATADIR}" --curve "${CURVE}" --output "${DATADIR}"/persistencebatch.csv + +# Calculate all the AoHs +littlejohn -j ${SLURM_JOB_CPUS_PER_NODE} -o "${DATADIR}"/aohbatch.log -c "${DATADIR}"/aohbatch.csv "${VIRTUAL_ENV}"/bin/python3 -- ./aoh-calculator/aohcalc.py --force-habitat + +# Generate validation summaries +python3 ./aoh-calculator/validation/collate_data.py --aoh_results "${DATADIR}"/aohs/current/ --output "${DATADIR}"/aohs/current.csv +python3 ./aoh-calculator/validation/collate_data.py --aoh_results "${DATADIR}"/aohs/pnv/ --output "${DATADIR}"/aohs/pnv.csv +for SCENARIO in "${SCENARIOS[@]}" +do + python3 ./aoh-calculator/validation/collate_data.py --aoh_results "${DATADIR}"/aohs/"${SCENARIO}"/ --output "${DATADIR}"/aohs/"${SCENARIO}".csv +done + +# Calculate predictors from AoHs +python3 ./aoh-calculator/summaries/species_richness.py --aohs_folder "${DATADIR}"/aohs/current/ \ + --output "${DATADIR}"/predictors/species_richness.tif +python3 ./aoh-calculator/summaries/endemism.py --aohs_folder "${DATADIR}"/aohs/current/ \ + --species_richness "${DATADIR}"/predictors/species_richness.tif \ + --output "${DATADIR}"/predictors/endemism.tif + +# Calculate the per species Delta P values +littlejohn -j ${SLURM_JOB_CPUS_PER_NODE} -o "${DATADIR}"/persistencebatch.log -c "${DATADIR}"/persistencebatch.csv "${VIRTUAL_ENV}"/bin/python3 -- ./deltap/global_code_residents_pixel.py + +for SCENARIO in "${SCENARIOS[@]}" +do + for TAXA in "${TAXAS[@]}" + do + python3 ./utils/raster_sum.py --rasters_directory "${DATADIR}"/deltap/"${SCENARIO}"/"${CURVE}"/"${TAXA}"/ --output "${DATADIR}"/deltap_sum/"${SCENARIO}"/"${CURVE}"/"${TAXA}".tif + done + + python3 ./utils/species_totals.py --aohs "${DATADIR}"/deltap/"${SCENARIO}"/"${CURVE}"/ --output "${DATADIR}"/deltap/"${SCENARIO}"/"${CURVE}"/totals.csv + + # Generate final map + python3 ./deltap/delta_p_scaled.py --input "${DATADIR}"/deltap_sum/"${SCENARIO}"/"${CURVE}"/ \ + --diffmap "${DATADIR}"/habitat/"${SCENARIO}"_diff_area.tif \ + --totals "${DATADIR}"/deltap/"${SCENARIO}"/"${CURVE}"/totals.csv \ + --output "${DATADIR}"/deltap_final/scaled_"${SCENARIO}"_"${CURVE}".tif +done diff --git a/utils/persistencegenerator.py b/utils/persistencegenerator.py index 633e82a..81f323a 100644 --- a/utils/persistencegenerator.py +++ b/utils/persistencegenerator.py @@ -1,31 +1,34 @@ import argparse +import sys from pathlib import Path import pandas as pd def species_generator( data_dir: Path, + curve: str, output_csv_path: Path, ): species_info_dir = data_dir / "species-info" taxas = [x.name for x in species_info_dir.iterdir()] + if curve not in ["0.1", "0.25", "0.5", "1.0", "gompertz"]: + sys.exit(f'curve {curve} not in expected set of values: ["0.1", "0.25", "0.5", "1.0", "gompertz"]') + res = [] for taxa in taxas: taxa_path = species_info_dir / taxa / "current" speciess = list(taxa_path.glob("*.geojson")) for scenario in ['restore_all', 'urban', 'arable', 'pasture', 'restore', 'restore_agriculture']: for species in speciess: - # for curve in ["0.1", "0.25", "0.5", "1.0", "gompertz"]: - for curve in ["0.25"]: - res.append([ - species, - data_dir / "aohs" / "current" / taxa, - data_dir / "aohs" / scenario / taxa, - data_dir / "aohs" / "pnv" / taxa, - curve, - data_dir / "deltap" / scenario / curve / taxa, - ]) + res.append([ + species, + data_dir / "aohs" / "current" / taxa, + data_dir / "aohs" / scenario / taxa, + data_dir / "aohs" / "pnv" / taxa, + curve, + data_dir / "deltap" / scenario / curve / taxa, + ]) df = pd.DataFrame(res, columns=[ '--speciesdata', @@ -47,6 +50,14 @@ def main() -> None: required=True, dest="data_dir", ) + parser.add_argument( + '--curve', + type=str, + choices=["0.1", "0.25", "0.5", "1.0", "gompertz"], + help='extinction curve, should be one of ["0.1", "0.25", "0.5", "1.0", "gompertz"]', + required=True, + dest="curve", + ) parser.add_argument( '--output', type=Path, @@ -56,7 +67,7 @@ def main() -> None: ) args = parser.parse_args() - species_generator(args.data_dir, args.output) + species_generator(args.data_dir, args.curve, args.output) if __name__ == "__main__": main()