Welcome to DOMINO-SEE¶
DOMINO-SEE is a framework for Detection Of Multi-layer Interconnected Network Occurrences for Spatial Extreme Events, leveraging event coincidences analysis and complex networks. It is built using xarray and can seamlessly benefit from the parallelization handling provided by dask. Its objective is to make it as simple as possible for users to construct event-based climate network analysis workflows.
The name DOMINO-SEE represents our approach to detecting and analyzing interconnected occurrences of hydroclimatic extreme events across spatial locations, inspired by the cascade effect of DOMINOes falling in a chain reaction. The SEE highlights the framework’s ability to capture the spatial synchronization and propagation of extreme events, emphasizing the interconnectedness inherent in complex environmental systems.
Key Features¶
Complex Network Generation: Fast and memory-efficient functions to build spatial networks from event series among spatial locations and multiple types (layers) of climate extreme events
Multi-dimensional Support: Native support for
xarray
DataArrays to handle multi-dimensional gridded climate dataParallel Processing:
dask
integration for efficient processing of large-scale climate datasetsBlockwise Computation: Utilities for splitting large spatial datasets into manageable blocks of netCDF datasets
Getting Started¶
To get started with dominosee, check out the Installation and Quickstart guides.
import numpy as np
import xarray as xr
# Create synthetic climate data
time = xr.cftime_range(start='2000-01-01', periods=365, freq='D')
lat = np.linspace(-90, 90, 10)
lon = np.linspace(-180, 180, 10)
# Generate random SPI-like data
np.random.seed(42)
spi = xr.DataArray(
np.random.normal(0, 1, (len(time), len(lat), len(lon))),
dims=('time', 'lat', 'lon'),
coords={'time': time, 'lat': lat, 'lon': lon},
name='SPI1'
)
# 1. Event detection
from dominosee.eventorize import get_event
da_event = get_event(spi, threshold=-1.0, extreme="below", event_name="drought")
# 2a. ECA analysis
from dominosee.eca import get_eca_precursor_from_events
da_precursor = get_eca_precursor_from_events(
eventA=da_event, eventB=da_event, delt=2, sym=True, tau=0
)
# 2b. ES analysis
from dominosee.eventorize import get_event_positions
from dominosee.es import get_event_sync_from_positions, create_null_model_from_indices
ds_event_pos = get_event_positions(da_event)
da_es = get_event_sync_from_positions(
positionsA=ds_event_pos.event_positions,
positionsB=ds_event_pos.event_positions,
time_dim='time',
tau_max=10
)
da_es_critical_value = create_null_model_from_indices(da_event.time, tau_max=10,
max_events=ds_event_pos.event_count.max(),
max_tau=10)
da_es_null = convert_null_model_for_locations(da_es_critical_value,
ds_event_pos.event_count,
ds_event_pos.event_count)
# 3. Network construction
from dominosee.network import get_link_from_confidence, get_link_from_critical_value
# For ECA
da_precursor_conf = get_eca_precursor_confidence(
precursor=da_precursor,
eventA=da_event,
eventB=da_event
)
da_precursor_link = get_link_from_confidence(da_precursor_conf, confidence_level=0.95)
# For ES (using a critical value threshold)
da_es_network = get_link_from_critical_values(da_es, critical_value=da_es_null, rule="greater_equal")
Examples¶
Check out the examples to see dominosee in action:
Event Coincidence Analysis (ECA) networks: Examples of using dominosee for ECA network analysis
Event Synchronization (ES) networks: Examples of using dominosee for ES network analysis
API Reference¶
The API Reference provides detailed documentation for all functions, classes, and methods in dominosee.
API Reference: Complete API documentation
Event Selection (dominosee.eventorize): Event selection
Event Coincidence Analysis (dominosee.eca): Event Coincidence Analysis (ECA)
Event Synchronization (dominosee.es): Event Synchronization (ES)
Network Construction (dominosee.network): Network creation