Computer-programs for the construction of seismicity models and maps based on Voronoi diagrams.
Warning
Note
For citation and a detailed presentation of the method, use the reference below:
Daniel, G., and Arroucau, P., 2025, Data-driven seismicity models based on Voronoi diagrams, Seismological Research Letters, vol. 96, Number 5, Sept. 2025, doi: 10.1785/0220240428.
Caution
Use this program at your own risks! We cannot guarantee the exacteness and the reliability of any program included in this repository.
This package consists in two main modules added of a suite of utilities for plotting or validation purposes:\
-
voronoi2density.py: Calculation of earthquake count/density grids based on Voronoi diagrams of a set of epicentral locations. This program can also propagate earthquake location and magnitude uncertainties using a random Monte-Carlo sampling process. This tends to produce a data-driven smoothing of spatial seismicity patterns. One count/density grid is constructed for each magnitude bin and is then projected onto a regular spatial mesh. -
compute_ab_values.py: Modelling of Frequency-Magnitude Distributions (FMD). This program uses a collection of count/density grids over several magnitude bins to reconstruct a FMD for each pixel of the regular spatial mesh. The program also estimates and returns a grid of parameters a and b from the Gutenberg & Richter law1 in each pixel. -
Plotting utilities:
plot_fmd.py: Plot the Frequency-Magnitude Distribution for one, or several, pixelsplot_map.py: Plot a map of earthquake counts/densities or a/b-values of the Gutenberg & Richter relationshipplot_distrib.py: Plot for one pixel of the regular grid the distribution of any quantity stored in GMT-formatted multiple segment files (i.e., counts, densities, a or b). This can be useful to inspect the distributions produced by the propagation of uncertainties usingvoronoi2density.py.
-
Testing utilities:
gof_tests.py: [ Under development ]. Goodness-of-fit tests based on residual analysis.- TODO : implement CSEP testing methods based on PyCSEP functionalites
TODO: Add sketches for the three module voronoi2density.py, compute_ab_values.py and make_plots.py showing:\
- Input required by each module
- Optional inputs
- Key role of the module
- Outputs produced by the module
- Optional outputs requested on-demand by the user
Notice for the propagation of magnitude uncertainties: When the user activates the propagation of locations & magnitude uncertainties, it is strongly advised to include events in the catalogue with magnitude smaller than the lower bound of the first bin in magnitude. This is important to avoid an exagerated depletion of earthquake counts in the first bin after perturbation of original magnitude values. Note however that, in any case, earthquake counts in each bin are expected to decrease due to the symmetrical (i.e., normal) perturbation of an exponentially-distributed variable (e.g., Rhoades, Tectonophysics, 1996).
Tutorials and examples will be included to assist interested users with our input formats, and with the use of the programs.
Syntax: python voronoi2density.py [options] parameters.txt
Note
A full list of available options can be obtained using the -h flag, i.e., python voronoi2density.py -h.
Estimation of (a, b) parameters of the Gutenberg & Richter law1
Syntax: python compute_ab_values.py [options] parameters.txt
Note
A full list of available options can be obtained using the -h flag, i.e., python compute_ab_values.py -h.
-
Configuration file (e.g.,
parameters.txt): Define the list of input files and settings for each application. Each line starting with "#" is considered as a comment and skipped.# VSM CONFIGURATION FILE # -- Define input files and output directories -- file_for_epicenters: epicenters_synthetics.txt file_for_geographical_bounds: bounds.txt file_for_magnitude_bins: bins.txt #file_for_FMD_limits_and_durations: fmd_info.txt #file_for_prior_b_information: b_prior_info.txt output_directory_for_files: results output_directory_for_figures: figures # -- Specification of input CRS and internal CRS used for the computation of areas -- input_CRS: EPSG:4326 # NB: EPSG:4326 --> WGS-84 internal_equal_area_CRS: EPSG:2154 # NB: EPSG:2154 --> Lambert 93 unit_for_internal_CRS_coordinates: m # -- Calculation settings -- mesh_discretization_step: 0.5 deg # NB: available units for step: km, deg density_scaling_factor: 1000.0 skip_ab_if_missing_priors: True define_completeness_automatically: False enable_verbosity: True # -- Options governing the propagation of uncertainties -- nb_bootstrap_samples: 100 perturb_magnitudes: True save_bootstrap_realizations: False nb_parallel_tasks: 3 b_value_to_remove_bias_on_perturbed_magnitudes: 1.0 save_bootstrap_realizations: False -
Magnitude bin information file (e.g.,
bins.txt): Definition of magnitude bins.
One bin per line, with columns arranged in the following order:[ID] [MIN] [MAX] [TMIN] [TMAX]. Each line starting with "#" is considered as a comment and skipped.# VSM MAGNITUDE BIN CONFIGURATION # [ID] [MIN] [MAX] [TMIN] [TMAX] 1 2.0 2.5 0.0 10000.0 2 2.5 3.0 0.0 10000.0 3 3.0 3.5 0.0 10000.0 4 3.5 4.0 0.0 10000.0 5 4.0 4.5 0.0 10000.0 6 4.5 5.0 0.0 10000.0 7 5.0 5.5 0.0 10000.0 8 5.5 6.0 0.0 10000.0 9 6.0 6.5 0.0 10000.0 10 6.5 7.0 0.0 10000.0 -
Target geographical area (e.g.,
bounds.txt): Define the target geographical area.
Bounds can be specified in two forms:- rectangular area, or
- any area enclosed in a specified polygon
Each line starting with "#" is considered as a comment and skipped.
Example format for a rectangular area: It is defined by its 4 corners, with columns arranged in
[LON] [LAT]order.# VSM TARGET GEOGRAPHICAL AREA -7.25 41.75 11.25 41.75 11.25 51.25 -7.25 51.25Example format for an area enclosed in a polygon:
It is defined by its 4 corners, with columns arranged in[LON] [LAT]order. The example below represents the polygon enclosing the Californian area used in the RELM2 earthquake forecasting experiment.# VSM TARGET GEOGRAPHICAL AREA -125.4 43.0 -118.16447368421052 43.0 -113.4 35.9 -113.1 35.45 -113.1 32.18333333333337 -113.2 31.8 -113.3 31.7 -113.5 31.6 -113.7 31.5 -118.36842105263159 31.5 -121.40 33.3 -121.5 33.4 -121.7 33.6 -122.0 33.9 -122.2 34.2 -123.2 35.7 -124.7 38.0 -125.4 39.1 -125.4 39.1 -125.4 43.0 -
Earthquake catalogue (e.g.,
epicenters.txt): Catalogue of earthquake epicentral locations with location and magnitude uncertainties.
One epicenter per line. Each line starting with "#" is considered as a comment and skipped.
The full format contains 9 columns, organized in the following order:[FLOATING DATE] [LON] [LAT] [MAG] [HALF-LENGTH OF SEMI-MAJOR AXIS IN KM] [HALF-LENGTH OF SEMI-MINOR AXIS IN KM] [AZIMUTH OF SEMI-MAJOR AXIS in DEG] [MAG UNCERTAINTY] [EVENT WEIGHT]# VSM EARTHQUAKE EPICENTERS # 1 EVENT PER LINE # COORDINATES EXPRESSED IN GEOGRAPHICAL COORDINATES (EPSG:4326) # LINE FORMAT: Floating_Date Longitude Latitude Magnitude Loc_Unc_SMAJ Loc_Unc_SMIN Loc_Unc_Az Mag_Unc Weight 1522.479452 0.75 46.917 5.8 100.0 100.0 0.0 0.43 1.0 1579.068493 2.0 46.583 6.0 50.0 50.0 0.0 0.47 1.0 1580.262295 1.5 51.0 5.8 50.0 50.0 0.0 0.49 1.0 1618.501370 -0.617 43.2 5.1 50.0 50.0 0.0 0.48 1.0 1640.510929 -1.367 48.933 5.2 100.0 100.0 0.0 0.42 1.0 1650.720548 7.6 47.55 5.3 50.0 50.0 0.0 0.44 1.0 1657.123288 0.617 47.117 5.7 50.0 50.0 0.0 0.51 1.0 1663.035616 -0.75 46.95 5.3 100.0 100.0 0.0 0.5 1.0 1665.093151 -0.05 43.1 5.4 100.0 100.0 0.0 0.43 1.0 1672.945355 7.75 47.483 5.6 50.0 50.0 0.0 0.45 1.0 1678.668493 5.783 43.75 5.0 100.0 100.0 0.0 0.4 1.0NB: Three alternative formats are available when uncertainties, or weights, are not specified in the catalogue:
- no uncertainties, no weights, 4-columns format, arranged in the following order:
[FLOATING DATE] [LON] [LAT] [MAG]. This implicitely sets all weights to 1.0. - with weights, but no uncertainties, 5-columns format, arranged in the following order:
[FLOATING DATE] [LON] [LAT] [MAG] [EVENT WEIGHT] - with uncertainties, but no weights, 8-columns format, arranged in the following order:
[FLOATING DATE] [LON] [LAT] [MAG] [HALF-LENGTH OF SEMI-MAJOR AXIS IN KM] [HALF-LENGTH OF SEMI-MINOR AXIS IN KM] [AZIMUTH OF SEMI-MAJOR AXIS in DEG] [MAG UNCERTAINTY]. This implicitely sets all weights to 1.0.
- no uncertainties, no weights, 4-columns format, arranged in the following order:
-
[Optional] FMD properties for each pixel (e.g.,
fmd_info.txt): Describes the bounds and completeness periods of frequency-magnitude distributions (FMD) for each pixel defined in the grid output by the programvoronoi2density.py. When completeness durations are specified in this file, they replace durations provided in the Magnitude Bin Configuration file, see above. One line per pixel, with semicolumn-delimited columns. Each line starting with "#" is considered as a comment and skipped.
Columns order: [CENTRAL LON] [CENTRAL LAT] [MMIN] [MMAX] [COMPLETENESS DURATION BIN 1] [COMPLETENESS DURATION BIN 2] ... [COMPLETENESS DURATION BIN N]
The first 3 columns are mandatory, others are optional.# VSM FMD INFORMATION FILE # LINE FORMAT: lon; lat; mmin; mmax; dur_1; dur_2; dur_3; dur_4; dur_5; dur_6; dur_7 -4.95; 43.85; 3.0; 6.5; 55.0; 60.0; 70.0; 220.0; 220.0; 220.0; 220.0 -4.95; 43.95; 3.0; 6.5; 55.0; 60.0; 70.0; 220.0; 220.0; 220.0; 220.0 -4.95; 44.05; 3.0; 6.5; 55.0; 60.0; 70.0; 220.0; 220.0; 220.0; 220.0 -4.95; 44.15; 3.0; 6.5; 55.0; 60.0; 70.0; 220.0; 220.0; 220.0; 220.0 -4.95; 44.25; 3.0; 6.5; 55.0; 60.0; 70.0; 220.0; 220.0; 220.0; 220.0Special values:
mmax = -9.0: Use an untruncated Gutenberg-Richter relationship for the pixel (i.e., equivalent tommax = inf)mmax = NaN: Use an untruncated Gutenberg-Richter relationship for the pixel (i.e., equivalent tommax = inf)- missing
mmax: Use an untruncated Gutenberg-Richter relationship for the pixel (i.e., equivalent tommax = inf)
-
[Optional] Prior b-value information for each pixel (e.g.,
b_prior_info.txt): Provides optional constraints on b-values for each pixel defined in the grid output by the programvoronoi2density.py. These constraints are expressed in terms of an a-priori normal distribution on b-values, parameterized by its mean and its standard-deviation. The estimation of (a, b) parameters for all pixels is realized by execution of the programcompute_ab_values.py, which loads this a-priori information only if the filename is listed in the configuration file. One line per pixel, with semicolumn-delimited columns. Each line starting with "#" is considered as a comment and skipped.
Columns order: [CENTRAL LON] [CENTRAL LAT] [B_MEAN] [B_STD_DEV]# B-PRIOR INFORMATION FILE # LINE FORMAT: lon; lat; b_mean; b_std -4.95; 47.75; 1.01; 0.1 -4.95; 47.85; 1.01; 0.1 -4.95; 47.95; 1.01; 0.1 -4.95; 48.05; 0.96; 0.1 -4.95; 48.15; 0.96; 0.1 -4.95; 48.25; 0.96; 0.1 -4.95; 48.35; 0.96; 0.1 -4.95; 48.45; 0.96; 0.1 -4.95; 48.55; 0.96; 0.1 -4.95; 48.65; 0.96; 0.1 -4.95; 48.75; 0.96; 0.1Special values:
b_mean = 0.0andb_std = 0.0: Flags indicating the absence of a-priori distribution for the pixelb_mean = -9.0andb_std = -9.0: Flags indicating the absence of a-priori distribution for the pixel
-
Earthquake counts per pixel (for each magnitude bin) (e.g.,
counts_bin_i.txt): GMT-formatted ASCII polygon file listing, for each pixel, the number of earthquake counts with magnitudes included in bin i. Each pixel is described by its 4 vertices. Z-values correspond to the estimated earthquake count within each pixel. -
Summary table of earthquake counts (for all bins and all pixels) (e.g.,
gridded_counts.txt): Table of earthquake count values for every bin and for every pixel in the target area, in CSV formatted, delimited by semi-colons (";"). One pixel per line. Number of columns is equal to the number of magnitude bins + 2.
Columns order: [CENTRAL LON] [CENTRAL LAT] [COUNTS IN BIN 1] [COUNTS IN BIN 2] ... [COUNTS IN BIN i] -
Earthquake density per pixel (for each magnitude bin) (e.g.,
density_bin_i.txt): GMT-formatted ASCII polygon file listing, for each pixel, the spatial density of earthquakes with magnitudes included in bin i. Each pixel is described by its 4 vertices. Density is obtained by dividing earthquake counts by pixel area. Z-values correspond to the spatial density of earthquake within each pixel. -
Summary table of earthquake densities (for all bins and all pixels) (e.g.,
gridded_densities.txt): Table of earthquake densities for every bin and for every pixel in the target area, in CSV formatted, delimited by semi-colons (";"). One pixel per line. Number of columns is equal to the number of magnitude bins + 2.
Columns order: [CENTRAL LON] [CENTRAL LAT] [DENSITY IN BIN 1] [DENSITY IN BIN 2] ... [DENSITY IN BIN i] -
Voronoi polygons (for each magnitude bin) (e.g.,
polygons_bin_i.txt): GMT-formatted ASCII polygon file listing all Voronoi polygons of the diagram obtained from epicentral locations of earthquakes with magnitudes included in bin i. Z-values correspond to the spatial density of earthquake within each polygon. This file is only produced when the random sampling of catalogue uncertainties is deactivated (i.e., optionnb_bootstrap_samplesset to 0). When uncertainties are propagated, a similar output can be produced for each random realisation by activating optionsave_bootstrap_realizations: Truein the Configuration file. -
Supplementary outputs when uncertainties are propagated using Monte-Carlo random sampling: The following files report standard deviations of earthquake count and density estimates when the option
nb_bootstrap_samples:is set with any number greater than 0 in the Configuration file.-
Standard-deviation of earthquake counts per pixel (for each magnitude bin) (e.g.,
counts_std_bin_i.txt): GMT-formatted ASCII polygon file listing, for each pixel, the standard-deviation of earthquake counts for events with magnitudes included in bin i. Each pixel is described by its 4 vertices. These files are only created when uncertainties are propagated by random sampling of catalogue uncertainties (i.e., when optionnb_bootstrap_samplesset to a number greater than 0). -
Standard-deviation of earthquake density per pixel (for each magnitude bin) (e.g.,
density_std_bin_i.txt): GMT-formatted ASCII polygon file listing, for each pixel, the standard-deviation of the spatial density of earthquakes with magnitudes included in bin i. Each pixel is described by its 4 vertices. These files are only created when uncertainties are propagated by random sampling of catalogue uncertainties (i.e., when optionnb_bootstrap_samplesset to a number greater than 0). -
Summary table of pixel-wise earthquake count standard deviations (for all bins) (e.g.,
gridded_counts_std.txt): Standard deviations of earthquake counts for every bin in each pixel. CSV format, same organisation than in filegridded_counts.txt -
Summary table of pixel-wise earthquake density standard deviations (for all bins) (e.g.,
gridded_densities_std.txt): Standard deviations of earthquake density for every bin in each pixel. CSV format, same organisation than in filegridded_densities.txt -
Optional outputs for each realization of the Monte-Carlo random sampling process: The next files store intermediary results produced during the Monte-Carlo random sampling process for the propagation of uncertainties.
For each realization, a new catalogue is built from the original catalogue, where earthquake locations (and optionally, magnitudes) are randomly drawn within the uncertainty domain specified in the original catalogue. Then, the process produces a Voronoi diagram based on this new catalogue, and it estimates earthquake counts and densities over the target gridded area from this diagram. The last step consists in averaging counts and densities over all realizations, and in storing final results in filesgridded_counts.txtandgridded_densities.txt.
Storing of intermediary results (catalogue, Voronoi polygons and pixel-wise estimates of earthquake counts and densities) can be activated by setting optionsave_bootstrap_realizations: True. By default, this option is deactivated (set toFalse).-
Random catalogue realization (e.g.,
bootstrap/catalog_bin_i_bs_j.txt): Random realization of the earthquake catalogue for bin i and for realization index j. This data is stored using a semi-colon (";")-delimited ASCII format with 4 columns.
Columns order: [TIME] [LONGITUDE] [LATITUDE] [MAGNITUDE] -
Earthquake counts per pixel (for each magnitude bin) for each realization (e.g.,
bootstrap/counts_bin_i_bs_j.txt): GMT-formatted ASCII polygon file listing, for each pixel, the number of earthquake counts with magnitudes included in bin i for realization index j. Same format than filecounts_bin_i.txt. -
Earthquake density per pixel (for each magnitude bin) for each realization (e.g.,
bootstrap/density_bin_i_bs_j.txt): GMT-formatted ASCII polygon file listing, for each pixel, the spatial density of earthquakes with magnitudes included in bin i for realization index j. Same format than filedensity_bin_i_bs_j.txt. -
Voronoi polygons (for each magnitude bin) for each realization (e.g.,
bootstrap/polygons_bin_i_bs_j.txt): GMT-formatted ASCII polygon file listing all Voronoi polygons obtained from epicentral locations of earthquakes in the random catalogue with magnitudes included in bin i for realization j. Same format than filepolygons_bin_i.txt.
-
-
[TO BE COMPLETED...]
[TO BE COMPLETED...]
Footnotes
-
Gutenberg, B., and Richter, C. F., 1944, Frequency of Earthquakes in California, Bulletin of the Seismological Society of America, 34, 4, pp.185-188. ↩ ↩2
-
Field, E. H., 2007, Overview of the Working Group for the Development of Regional Earthquake Likelihood Models (RELM), Seismological Research Letters, 78, 1, pp.7-16 ↩