# Important: If you haven't already, you need to install pygbif using the command below
# Run this command in the terminal to install: pip install pygbif

# Import packages
import os # For file paths, data management, interacting with the operating system
import pathlib # Object oriented file path management
import time # Time how long operations take in Python
import zipfile # Deal with zip files
from getpass import getpass # Security login for downloading GBIF data
from glob import glob # Search for specific file extensions

# Libraries for data analysis 
import pandas as pd # Data analysis for tabular data
import geopandas as gpd # Data analysis for geospatial data
import pygbif.occurrences as occ # Get occurrence data
import pygbif.species as species # Get species data
import requests # Download data

# Libraries for Dynamic mapping
import cartopy.crs as ccrs #Coordinate reference systems 
import panel as pn #customize the dynamic map
import hvplot.pandas #make the holoview plot with geo-dataframe
import calendar # Change month numbers to month names

# Define data directory in user's home directory
data_dir_2004 = os.path.join(
    # Home directory
    pathlib.Path.home(),
    # Earth analytics data directory
    'earth-analytics',
    'data',
    # Project directory
    'redstart-migration-2004',
)
print(data_dir_2004)
# Make the directory using the path defined above
os.makedirs(data_dir_2004, exist_ok=True)

C:\Users\warno\earth-analytics\data\redstart-migration-2004

# GBIF LOGIN CODE CHUNK
# Code credit: Elsa Culler
# This code allows us to download data in Python from gbif.org
# This code ASKS for your credentials for gbif.org
# and saves it for the rest of the session.
# NEVER put your credentials into your code

# GBIF needs a username, password, and email 
# All 3 need to match the account
reset = False

# Request and store username
if (not ('GBIF_USER'  in os.environ)) or reset:
    os.environ['GBIF_USER'] = input('GBIF username:')

# Securely request and store password
if (not ('GBIF_PWD'  in os.environ)) or reset:
    os.environ['GBIF_PWD'] = getpass('GBIF password:')
    
# Request and store account email address
if (not ('GBIF_EMAIL'  in os.environ)) or reset:
    os.environ['GBIF_EMAIL'] = input('GBIF email:')

# Search for the American Redstart (Setophaga ruticilla)
backbone = species.name_backbone(name='Setophaga ruticilla')

# Save the species key
species_key = backbone["speciesKey"]
backbone

{'usageKey': 2489985,
 'scientificName': 'Setophaga ruticilla (Linnaeus, 1758)',
 'canonicalName': 'Setophaga ruticilla',
 'rank': 'SPECIES',
 'status': 'ACCEPTED',
 'confidence': 99,
 'matchType': 'EXACT',
 'kingdom': 'Animalia',
 'phylum': 'Chordata',
 'order': 'Passeriformes',
 'family': 'Parulidae',
 'genus': 'Setophaga',
 'species': 'Setophaga ruticilla',
 'kingdomKey': 1,
 'phylumKey': 44,
 'classKey': 212,
 'orderKey': 729,
 'familyKey': 5263,
 'genusKey': 2489984,
 'speciesKey': 2489985,
 'class': 'Aves'}

gbif_pattern_2004 = os.path.join(data_dir_2004, '*.csv')
if not glob(gbif_pattern_2004):
    # Only submit one request
    if not 'GBIF_DOWNLOAD_KEY' in os.environ:
        # Submit query to GBIF
        gbif_query = occ.download([
            f"speciesKey = { species_key }",
            "hasCoordinate = True",
            "year = 2004",
        ])
        os.environ['GBIF_DOWNLOAD_KEY'] = gbif_query[0]
        print("Submitted query to GBIF")
    
    # Wait for download to build
    download_key = os.environ['GBIF_DOWNLOAD_KEY']
    wait = occ.download_meta(download_key)['status']
    while not  wait=='SUCCEEDED':
        wait = occ.download_meta(download_key)['status']
        time.sleep(5)
    print("Building download")

    # Download GBIF data
    print("Downloading GBIF Data")
    download_info = occ.download_get(
        os.environ['GBIF_DOWNLOAD_KEY'],
        path=data_dir_2004)
    print("Finished Downloading GBIF Data")
    
    # Unzip GBIF data
    with zipfile.ZipFile(download_info['path']) as download_zip:
        download_zip.extractall(path=data_dir_2004)
    print("Unzipping GBIF data")

# Find extracted .csv file path and take 1st result
gbif_path_2004 = glob(gbif_pattern_2004)[0]

with open(gbif_path_2004) as f:
    for i in range(2):
        print(f.readline().strip())

gbifID	datasetKey	occurrenceID	kingdom	phylum	class	order	family	genus	species	infraspecificEpithet	taxonRank	scientificName	verbatimScientificName	verbatimScientificNameAuthorship	countryCode	locality	stateProvince	occurrenceStatus	individualCount	publishingOrgKey	decimalLatitude	decimalLongitude	coordinateUncertaintyInMeters	coordinatePrecision	elevation	elevationAccuracy	depth	depthAccuracy	eventDate	day	month	year	taxonKey	speciesKey	basisOfRecord	institutionCode	collectionCode	catalogNumber	recordNumber	identifiedBy	dateIdentified	license	rightsHolder	recordedBy	typeStatus	establishmentMeans	lastInterpreted	mediaType	issue
985926802	4fa7b334-ce0d-4e88-aaae-2e0c138d049e	URN:catalog:CLO:EBIRD:OBS134172612	Animalia	Chordata	Aves	Passeriformes	Parulidae	Setophaga	Setophaga ruticilla		SPECIES	Setophaga ruticilla (Linnaeus, 1758)	Setophaga ruticilla		US	Rose City	Michigan	PRESENT	5	e2e717bf-551a-4917-bdc9-4fa0f342c530	44.42142	-84.1167							2004-07-25	25	7	2004	2489985	2489985	HUMAN_OBSERVATION	CLO	EBIRD	OBS134172612				CC_BY_4_0		obsr31360			2025-11-05T07:38:56.818Z		CONTINENT_DERIVED_FROM_COORDINATES;TAXON_CONCEPT_ID_NOT_FOUND

# Load the GBIF data, this time defining the delimiter (tabular) represents as '\t'
# We also set the index to gbifID, and select the columns of interest (lat, long, and month)
gbif_df = pd.read_csv(
    gbif_path_2004,
    delimiter='\t',
    index_col='gbifID',
    usecols=['gbifID','decimalLongitude','decimalLatitude','month']
)
gbif_df.head()

# Converting the dataframe to a geodataframe
# Using geopandas library
gbif_gdf = (
    gpd.GeoDataFrame(
        gbif_df, 
        geometry=gpd.points_from_xy(
            gbif_df.decimalLongitude, 
            gbif_df.decimalLatitude), 
        crs="EPSG:4326") # Manually set the crs to match the ecoregions data
    # Select the desired columns
    [['month','geometry']]
)
gbif_gdf

# Ecoregions data from GeographyRealm, RESOLVE, and OneEarth: 
# GeographyRealm: https://www.geographyrealm.com/terrestrial-ecoregions-gis-data/
# OneEarth: https://www.oneearth.org/announcing-the-release-of-ecoregion-snapshots/
eco_url = (
    "https://storage.googleapis.com/teow2016/Ecoregions2017.zip")

# Set up a path to save the data to your machine
eco_dir = os.path.join(data_dir_2004, 'resolve_ecoregions') 
# Make ecoregions directory
os.makedirs(eco_dir, exist_ok=True)

# Join ecoregions shapefile path
eco_path = os.path.join(eco_dir, 'ecoregions.shp')

# Only download once
if not os.path.exists(eco_path):
    my_gdf = gpd.read_file(eco_url)
    my_gdf.to_file(eco_path)

# Opening the download and renaming some of the column names
eco_gdf = (
    gpd.read_file(eco_path)
    [['OBJECTID','ECO_NAME','SHAPE_AREA','geometry']]
    .rename(columns={
        'OBJECTID': 'ecoregion_id',
        'ECO_NAME': 'name',
        'SHAPE_AREA': 'area'
    })
    .set_index('ecoregion_id')
)

# Plot the ecoregions quickly to check download
eco_gdf.plot()
eco_gdf.head()

gbif_ecoregion_gdf = (
    eco_gdf
    # Match the CRS of the GBIF data and the ecoregions
    .to_crs(gbif_gdf.crs)
    # Find ecoregion for each observation
    .sjoin(
        gbif_gdf,
        how='inner', 
        predicate='contains')
)
gbif_ecoregion_gdf

# Count the observations in each ecoregion each month
occurrence_df = (
    gbif_ecoregion_gdf
    # Reset index
    .reset_index()
    # For each ecoregion, for each month...
    .groupby(['ecoregion_id','month'])
    # ...count the number of occurrences
    # creating a new column 'occurrences' that counts all the gbif IDs in a group
    .agg(
        occurrences=('gbifID', 'count'),
        area=('area','first')) 
)

#Normalize by area
occurrence_df['density'] = (
    occurrence_df.occurrences
    / occurrence_df.area
)

# Get rid of rare observations (possible misidentification?)
occurrence_df = occurrence_df[occurrence_df.occurrences>1]

# Take the mean by ecoregion
mean_occurrences_by_ecoregion = (
    occurrence_df
    .groupby('ecoregion_id')
    .mean()
)
mean_occurrences_by_ecoregion

# Take the mean by month
mean_occurrences_by_month = (
    occurrence_df
    .groupby('month')
    .mean()
)
mean_occurrences_by_month

# Normalize by space and time for sampling effort
occurrence_df['norm_occurrences'] = (
    occurrence_df[['density']]
    / mean_occurrences_by_ecoregion[['density']]
    / mean_occurrences_by_month[['density']]
)
occurrence_df

occurrence_df.norm_occurrences.plot.hist()

<Axes: ylabel='Frequency'>

# Plotting a scatter with regular "occurrences" by month
occurrence_df.reset_index().plot.scatter(
    x='month', y='occurrences', c='ecoregion_id', logy=True
)

<Axes: xlabel='month', ylabel='occurrences'>

# Plotting a scatter with normalized occurrences by month
occurrence_df.reset_index().plot.scatter(
    x='month', y='norm_occurrences', c='ecoregion_id', logy=True
)

<Axes: xlabel='month', ylabel='norm_occurrences'>

# Simplify the geometry to speed up processing
eco_gdf.geometry = eco_gdf.simplify(.1, preserve_topology=False)

# Change the CRS to Mercator for mapping
eco_gdf = eco_gdf.to_crs(ccrs.Mercator())

# Check that the plot runs in a reasonable amount of time
eco_gdf.hvplot(geo=True, crs=ccrs.Mercator())

# Join the occurrences with the plotting GeoDataFrame
occurrence_gdf = eco_gdf.join(occurrence_df[['norm_occurrences']])

# Remove any empty geometry rows. These can cause errors if they persist in the data. 
occurrence_gdf = occurrence_gdf[~occurrence_gdf.geometry.is_empty]

# Get the plot bounds so they don't change with the slider
# xmin, ymin, xmax, ymax = gbif_gdf.to_crs(ccrs.Mercator()).total_bounds

# This code above will plot the entire range of the data (even other continents/unusual sightings)
# I am just interested in the South-North American data, so I set the bounds manually.

xmin = -35000000.0
xmax = -2000000.0
ymax = 11000000.0
ymin = -4000000.0

# List month names on the slider
month_widget = pn.widgets.DiscreteSlider(
    options={
        calendar.month_name[month_num]: month_num
        for month_num in range(1,13)
    }
)

# Plot occurrence by ecoregion and month
migration_plot = (
    occurrence_gdf
    .hvplot(
        c='norm_occurrences',
        groupby='month',
        # Use background tiles
        geo=True, crs=ccrs.Mercator(), tiles='CartoLight',
        title="American Redstart Migration, 2004",
        xlim=(xmin, xmax), ylim=(ymin, ymax),
        frame_height=220,
        frame_width=220,
        widgets={'month': month_widget},
        widget_location='bottom'
    )
)

# Save the plot
migration_plot.save('redstart-migration-2004.html', embed=True)
# Display the plot
migration_plot

WARNING:bokeh.core.validation.check:W-1005 (FIXED_SIZING_MODE): 'fixed' sizing mode requires width and height to be set: figure(id='p11490', ...)

BokehModel(combine_events=True, render_bundle={'docs_json': {'9b0f8a9d-40cc-4b2c-b432-1d021b3acbac': {'version…

occurrence_gdf.to_file("redstart_occurrence_2004.geojson", driver="GeoJSON")

INFO:Created 440 records

	decimalLatitude	decimalLongitude	month
gbifID
985926802	44.421420	-84.116700	7.0
985836732	36.553740	-89.328010	10.0
985770190	36.577732	-88.147026	10.0
985746922	36.642014	-84.321660	9.0
985519421	22.326576	-81.188965	4.0

	month	geometry
gbifID
985926802	7.0	POINT (-84.1167 44.42142)
985836732	10.0	POINT (-89.32801 36.55374)
985770190	10.0	POINT (-88.14703 36.57773)
985746922	9.0	POINT (-84.32166 36.64201)
985519421	4.0	POINT (-81.18896 22.32658)
...	...	...
1147139514	5.0	POINT (-76.4532 42.48045)
1147116162	5.0	POINT (-74.55236 44.17583)
1147108424	6.0	POINT (-76.4532 42.48045)
1147107409	6.0	POINT (-74.66193 44.2033)
1147035591	5.0	POINT (-76.4532 42.48045)

	name	area	geometry
ecoregion_id
1.0	Adelie Land tundra	0.038948	MULTIPOLYGON (((158.7141 -69.60657, 158.71264 ...
2.0	Admiralty Islands lowland rain forests	0.170599	MULTIPOLYGON (((147.28819 -2.57589, 147.2715 -...
3.0	Aegean and Western Turkey sclerophyllous and m...	13.844952	MULTIPOLYGON (((26.88659 35.32161, 26.88297 35...
4.0	Afghan Mountains semi-desert	1.355536	MULTIPOLYGON (((65.48655 34.71401, 65.52872 34...
5.0	Ahklun and Kilbuck Upland Tundra	8.196573	MULTIPOLYGON (((-160.26404 58.64097, -160.2673...

	name	area	geometry	gbifID	month
ecoregion_id
13.0	Alberta-British Columbia foothills forests	17.133639	MULTIPOLYGON (((-119.53979 55.81661, -119.5443...	3214209708	6.0
13.0	Alberta-British Columbia foothills forests	17.133639	MULTIPOLYGON (((-119.53979 55.81661, -119.5443...	5674975578	6.0
13.0	Alberta-British Columbia foothills forests	17.133639	MULTIPOLYGON (((-119.53979 55.81661, -119.5443...	5523202791	6.0
13.0	Alberta-British Columbia foothills forests	17.133639	MULTIPOLYGON (((-119.53979 55.81661, -119.5443...	256786965	6.0
13.0	Alberta-British Columbia foothills forests	17.133639	MULTIPOLYGON (((-119.53979 55.81661, -119.5443...	5551786796	6.0
...	...	...	...	...	...
839.0	Northern Rockies conifer forests	35.905513	POLYGON ((-119.99977 54.53117, -119.8914 54.45...	1784040249	7.0
839.0	Northern Rockies conifer forests	35.905513	POLYGON ((-119.99977 54.53117, -119.8914 54.45...	138754461	6.0
839.0	Northern Rockies conifer forests	35.905513	POLYGON ((-119.99977 54.53117, -119.8914 54.45...	138753536	9.0
844.0	Lesser Antillean dry forests	0.053407	MULTIPOLYGON (((-61.65263 12.23089, -61.63773 ...	5606006408	5.0
844.0	Lesser Antillean dry forests	0.053407	MULTIPOLYGON (((-61.65263 12.23089, -61.63773 ...	5596170282	4.0

	occurrences	area	density
ecoregion_id
13.0	12.000000	17.133639	0.700377
17.0	60.800000	7.958751	7.639390
23.0	2.500000	3.346216	0.747113
33.0	74.714286	16.637804	4.490634
34.0	33.428571	18.674884	1.790028
...	...	...	...
827.0	4.000000	9.664680	0.413878
828.0	7.000000	64.674744	0.108234
833.0	2.000000	0.610793	3.274429
838.0	6.500000	4.286144	1.516515
839.0	19.200000	35.905513	0.534737

	occurrences	area	density
month
1.0	7.666667	5.683188	12.988499
2.0	7.771429	5.180762	51.662875
3.0	5.903226	6.152404	29.587264
4.0	17.628571	8.718374	19.557518
5.0	84.436364	18.996233	14.146127
6.0	68.500000	26.174810	3.453364
7.0	32.428571	29.147840	1.581037
8.0	26.578947	24.332232	2.349243
9.0	40.782609	19.425837	21.520425
10.0	18.725000	11.044585	49.997778
11.0	5.500000	7.004568	4.030330
12.0	6.363636	5.405988	20.823892

GBIF Data Download and Methods for 2004¶