{'fraction_non_alphanumeric': 0.025937760888965655, 'fraction_numerical': 0.009902985289460806, 'mean_word_length': 4.821602145131557, 'pattern_counts': {'":': 0, '<': 4, '<?xml version=': 0, '>': 5, 'https://': 1, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 1, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Earthquake magnitude and location estimation from real time seismic waveforms with a transformer network

Precise real time estimates of earthquake magnitude and location are essential for early warning and rapid response. While recently multiple deep learning approaches for fast assessment of earthquakes have been proposed, they usually rely on either seismic records from a single station or from a fixed set of seismic stations. Here we introduce a new model for real-time magnitude and location estimation using the attention based transformer networks. Our approach incorporates waveforms from a dynamically varying set of stations and outperforms deep learning baselines in both magnitude and location estimation performance. Furthermore, it outperforms a classical magnitude estimation algorithm considerably and shows promising performance in comparison to a classical localization algorithm. In this work, we furthermore conduct a comprehensive study of the requirements on training data, the training procedures and the typical failure modes using three diverse and large scale data sets. Our analysis gives several key insights. First, we can precisely pinpoint the effect of large training data; for example, a four times larger training set reduces the required time for real time assessment by a factor of four. Second, the basic model systematically underestimates large magnitude events. This issue can be mitigated by incorporating events from other regions into the training through transfer learning. Third, location estimation is highly precise in areas with sufficient training data, but is strongly degraded for events outside the training distribution. Our analysis suggests that these characteristics are not only present for our model, but for most deep learning models for fast assessment published so far. They result from the black box modeling and their mitigation will likely require imposing physics derived constraints on the neural network.


Introduction
Recently, multiple studies investigated deep learning on raw seismic waveforms for the fast assessment of earthquake parameters, such as magnitude (e.g. Lomax et al. (2019), Mousavi & Beroza (2020), van den Ende & Ampuero (2020)), location (e.g. Kriegerowski et al. (2019), Mousavi & Beroza (2019), van den Ende & Ampuero (2020)) and peak ground acceleration (e.g. Jozinović et al. (2020)). Deep learning is well suited for these tasks, as it does not rely on manually selected features, but can learn to extract relevant information from the raw input data. This property allows the models to use the full information contained in the waveforms of an event. However, the models published so far use fixed time windows and can not be applied to data of varying length without retraining. Similarly, except the model by van den arXiv:2101.02010v1 [physics.geo-ph] 6 Jan 2021 Ende & Ampuero (2020), all models process either waveforms from only a single seismic station or rely on a fixed set of seismic stations defined at training time. The model by van den Ende & Ampuero (2020) enables the use of a variable station set, but combines measurements from multiple stations using a simple pooling mechanism. While it has not been studied so far in a seismological context, it has been shown in the general domain that set pooling architectures are in practice limited in the complexity of functions they can model (Lee et al., 2019).
Here we introduce a new model for magnitude and location estimation based on the architecture recently introduced for the transformer earthquake alerting model (TEAM) (Münchmeyer et al., 2020a), a deep learning based earthquake early warning model. While TEAM estimated PGA at target locations, our model estimates magnitude and hypocentral location of the event. We call our adaptation TEAM-LM, TEAM for location and magnitude estimation. We use TEAM as a basis due to its flexible multi-station approach and its ability to process incoming data effectively in real-time, issuing updated estimates as additional data become available. Similar to TEAM, TEAM-LM uses mixture density networks to provided probability distributions rather than merely point estimates as predictions. For magnitude estimation, our model outperforms two state of the art baselines, one using deep learning (van den Ende & Ampuero, 2020) and one classical approach (Kuyuk & Allen, 2013). For location estimation, our model outperforms a deep learning baseline (van den Ende & Ampuero, 2020) and shows promising performance in comparison to a classical localization algorithm.
We note a further deficiency of previous studies for deep learning in seismology. Many of these pioneering studies focused their analysis on the average performance of the proposed models. Therefore, little is known about the conditions under which these models fail, the impact of training data characteristics, the possibility of sharing knowledge across world regions, and of specific training strategies. All of these are of particular interest when considering practical application of the models.
To address these issues and provide guidance for practitioners, we perform a comprehensive evaluation of TEAM-LM on three large and diverse data sets: a regional broadband data set from Northern Chile, a strong motion data set from Japan, and another strong motion data set from Italy. These data sets differ in their seismotectonic environment (North Chile and Japan: subduction zones; Italy: dominated by both convergent and divergent continental deformation), their spatial extent (North Chile: regional scale; Italy and Japan: national catalogs), and the instrument type (North Chile: broadband, Italy and Japan: strong motion). All three data sets contain hundreds of thousands of waveforms. North Chile is characterized by a relatively sparse station distribution, but a large number of events and a low magnitude of completeness. There are far more stations in the Italy and Japan data sets, but a smaller number of earthquakes. This selection of diverse data sets allows for a comprehensive analysis, giving insights for different use cases. Our targeted experiments show that the characteristics are rooted in the principle structure used by TEAM-LM, i.e., the black box approach of learning a very flexible model from data, without imposing any physical constraints. As this black box approach is common to all current fast assessment models using deep learning, they can be transferred to these models. This finding is further backed by comparison to the results from previous studies.

Datasets
For this study we use three data sets (Table 1, Figure 1): one from Northern Chile, one from Italy and one from Japan. The Chile data set is based on the catalog by Sippl et al. (2018) with the magnitude values from Münchmeyer et al. (2020b). While there were minor changes in the seismic network configuration during the time covered by the catalog, the station set used in the construction of this catalog had been selected to provide a high degree of stability of the locations accuracy throughout the observational period (Sippl et al., 2018). Similarly, the magnitude scale has been carefully calibrated to achieve a high degree of consistency in spite of significant variations of attenuation (Münchmeyer et al., 2020b). This data set therefore contains the highest quality labels among the data sets in this study. For the Chile data set, we use broadband seismogramms from the fixed set of 24 stations used for the creation of the original catalog and magnitude scale. Although the Chile data set has the smallest number of stations of the three data sets, it comprises three to four times as many waveforms as the other two due to the large number of events.
The data sets for Italy and Japan are more focused on early warning, containing fewer events and only strong motion waveforms. They are based on catalogs from the INGV (ISIDe Working Group, 2007) and the NIED KiKNet (National Research Institute For Earth Science And Disaster Resilience, 2019), respectively. The data sets each encompass a larger area than the Chile data set and include waveforms from significantly more stations. In contrast to the Chile data sets, the station coverage differs strongly between different events, as only stations recording the event are considered. In particular, KiKNet stations do not record continuous waveforms, but operate in trigger mode, only saving waveforms if an event triggered at the station. For Japan each station comprises two sensors, one at the surface and one borehole sensor. Therefore for Japan we have 6 component recordings (3 surface, 3 borehole) available instead of the 3 component  Table 1: Overview of the data sets. The lower boundary of the magnitude category is the 5th percentile of the magnitude; this limit is chosen as each data set contains a small number of unrepresentative very small events. The upper boundary is the maximum magnitude. Magnitudes are given with two digit precision for Chile, as the precision of the underlying catalog is higher than for Italy and Japan. The Italy data set uses different magnitudes for different events, which are M L (>90 % of the events), M W (<10 %) and m b (<1 %). For depth and distance minimum, median and maximum are stated. Distance refers to the epicentral distance between stations and events. Note that the count of traces refers to the number of waveform-triplets (for three components, or group of six waveforms for the Japanese stations). The sensor types are broadband (BB) and strong motion (SM).
For each data set we use the magnitude scale provided in the catalog. For the Chile catalog, this is M A , a peak displacement based scale, but without the Wood-Anderson response and therefore saturation-free for large events (Münchmeyer et al., 2020b;Deichmann, 2018). For Japan M JMA is used. M JMA combines different magnitude scales, but similarly to M A primarily uses horizontal peak displacement (Doi, 2014). For Italy the catalog provides different magnitude types approximately dependent on the size of the event: M L (>90 % of the events), M W (<10 %) and m b (<1 %). We note that while the primary magnitude scales for all data sets are peak-displacement based, the precision of the magnitudes vary, with the highest precision for Chile. This might lead to slightly worse magnitude estimation performance for Italy and Japan.
For all data sets the data were not subselected based on the type of seismicity but only based on the location (for Chile and Italy) or depending if they triggered (Japan). This guarantees that, even though we made use of a catalog to assemble our training data, the resulting data sets are suitable for training and assessing methods geared at real-time applications without any prior knowledge about the earthquakes. We focus on earthquake characterization and do not discuss event detection or separation from noise; we refer the interested reader to Perol et al. (2018); Mousavi et al. (2018).
We split each data set into training, development and test set. For Chile and Japan we apply a simple chronological split with approximate ratios of 60:10:30 between training, development and test set, with the most recent events in the test set. As the last 30% of the Italy data set consist of less interesting events for early warning, we instead use all events from 2016 as test set and the remaining events as training and development sets. We reserve all of 2016 for testing, as it contains a long seismic sequence in central Italy with two mainshocks in August (M W = 6.5) and October (M W = 6.0). Notably, the largest event in the test set is significantly larger than the largest event in the training set (M w = 6.1 L'Aquila event in 2007), representing a challenging test case. For Italy, we assign the remaining events to training and development set randomly with a 6:1 ratio.
The Italy data set has been published as Münchmeyer et al. (2020). Download instructions for the Japan data set are available at https://github.com/yetinam/TEAM. A data publication for the Chile data set with GFZ data services is under preparation.

The transformer earthquake alerting model for magnitude and location
We build a model for real time earthquake magnitude and location estimation based on the core ideas of the transformer earthquake alerting model (TEAM), as published in Münchmeyer et al. (2020a). TEAM is an end-to-end peak ground Figure 2: Overview of the adapted transformer earthquake alerting model, showing the input, the feature extraction, the feature combination, the magnitude/location estimation and the output. For simplicity, not all layers are shown, but only their order and combination is visualized schematically. For the exact number of layers and the size of each layer please refer to tables SM 2 to SM 4. Please note that the number of input stations is variable, due to the self-attention mechanism in the feature combination.
acceleration (PGA) model calculating probabilistic PGA estimates based on incoming waveforms from a flexible set of stations. It employs the transformer network method (Vaswani et al., 2017), an attention based neural network which was developed in the context of natural language processing (NLP), at the core of its algorithm. Here, we adapt TEAM to calculate real time probabilistic estimates of event magnitude and hypocentral location. As our model closely follows the architecture and key ideas of TEAM, we use the name TEAM-LM to refer to the location and magnitude estimation model.
Similar to TEAM, TEAM-LM consists of three major components ( Figure 2): a feature extraction, which generates features from raw waveforms at single stations, a feature combination, which aggregates features across multiple stations, and an output estimation. Here, we briefly discuss the core ideas of the TEAM architecture and training and put a further focus on the necessary changes for magnitude and location estimation. For a more detailed account of TEAM and TEAM-LM we refer to Münchmeyer et al. (2020a), Tables SM 2 to SM 4 and the published implementation.
The input to TEAM consists of three component seismogramms from multiple stations and their locations. TEAM aligns all seismogramms to start and end at the same times t 0 and t 1 . We choose t 0 to be 5 seconds before the first P arrival at any station. This allows the model to understand the noise conditions at all stations. We limit t 1 to be at latest t 0 + 30 s. In a real-time scenario t 1 is the current time, i.e., the available amount of waveforms, and we use the same approach to imitate real-time waveforms in training and evaluation. The waveforms are padded with zeros to a length of 30 s to achieve constant length input to the feature extraction.
TEAM uses a CNN architecture for feature extraction, which is applied separately at each station. The architecture consists of several convolution and pooling layers, followed by a multi-layer perceptron (Table SM 2). To avoid scaling issues, each input waveform is normalized through division by its peak amplitude. As the amplitude is expected to be a key predictor for the event magnitude, we provide the logarithm of the peak amplitude as a further input to the multi-layer perceptron inside the feature extraction network. We ensure that this transformation does not introduce a knowledge leak by calculating the peak amplitude only based on the waveforms until t 1 . The full feature extraction returns one vector for each station, representing the measurements at the station.
The feature vectors from multiple stations are combined using a transformer network (Vaswani et al., 2017). Transformers are attention based neural networks, originally introduced for natural language processing. A transformer takes a set of n vectors as input, and outputs again n vectors which now incorporate the context of each other. The attention mechanism allows the transformer to put special emphasis on inputs that it considers particularly relevant and thereby model complex inter-station dependencies. Importantly, the parameters of the transformer are independent of the number of input vectors n, allowing to train and apply a transformer on variable station sets. To give the transformer a notion of the position of the stations, TEAM encodes the latitude, longitude and elevation of the stations using a sinusoidal embedding and adds this embedding to the feature vectors.
TEAM adds the position embeddings of the PGA targets as additional inputs to the transformer. In TEAM-LM, we aim to extract information about the event itself, where we do not know the position in advance. To achieve this, we add an event token, which is a vector with the same dimensionality as the positional embedding of a station location, and which can be thought of as a query vector. This approach is inspired by the so-called sentence tokens in NLP that are used to extract holistic information on a sentence (Devlin et al., 2018). The elements of this event query vector are learned during the training procedure.
From the transformer output, we only use the output corresponding to the event token, which we term event embedding and which is passed through another multi-layer perceptron predicting the parameters and weights of a mixture of Gaussians (Bishop, 1994). We use N = 5 Gaussians for magnitude and N = 15 Gaussians for location estimation. For computational and stability reasons, we constrain the covariance matrix of the individual Gaussians for location estimation to a diagonal matrix to reduce the output dimensionality. Even though uncertainties in latitude, longitude and depth are known to generally be correlated, this correlation can be modeled with diagonal covariance matrices by using the mixture.
The model is trained end-to-end using a log-likelihood loss with the Adam optimizer (Kingma & Ba, 2014). We train separate models for magnitude and for location. As we observed difficulties in the onset of the optimization when starting from a fully random initialization, we pretrain the feature extraction network. To this end we add a mixture density network directly after the feature extraction and train the resulting network to predict magnitudes from single station waveforms. We then discard the mixture density network and use the weights of the feature extraction as initialization for the end-to-end training. We use this pretraining method for both magnitude and localization networks.
Similarly to the training procedure for TEAM we make extensive use of data augmentation during training. First, we randomly select a subset of up to 25 stations from the available station set. We limit the maximum number to 25 for computational reasons. Second, we apply temporal blinding, by zeroing waveforms after a random time t 1 . This type of augmentation allows TEAM-LM to be applied to real time data. To avoid knowledge leaks for Italy and Japan, we only use stations as inputs that triggered before time t 1 for these data sets. This is not necessary for Chile, as there the maximum number of stations per event is below 25 and waveforms for all events are available for all stations active at that time, irrespective of whether the station actually recorded the event. Third, we oversample large magnitude events, as they are strongly underrepresented in the training data set. We discuss the effect of this augmentation in further detail in the Results section. In contrast to the station selection during training, in evaluation we always use the 25 stations picking first. Again, we only use stations and their waveforms as input once they triggered, thereby ensuring that the station selection does not introduce a knowledge leak.

Baseline methods
Recently, van den Ende & Ampuero (2020) suggested a deep learning method capable of incorporating waveforms from a flexible set of stations. Their architecture uses a similar CNN based feature extraction as TEAM-LM. In contrast to TEAM-LM, for feature combination it uses maximum pooling to aggregate the feature vectors from all stations instead of a transformer. In addition they do not add predefined position embeddings, but concatenate the feature vector for each station with the location coordinates and apply a multi-layer perceptron to get the final feature vectors for each station. The model of van den Ende & Ampuero (2020) is both trained and evaluated on 100 s long waveforms. Therefore it is not suitable for real time processing. Furthermore, the authors did not publish an implementation of the model. For these reasons we are not able to compare the exact model of van den Ende & Ampuero (2020) to TEAM-LM.
To still compare TEAM-LM to the techniques introduced in this approach, we implemented a model based on the key concepts of van den Ende & Ampuero (2020). As we aim to evaluate the performance differences from the conceptual changes, rather than different hyperparameters, e.g., the exact size and number of the convolutional layers, we use the same architecture as TEAM-LM for the feature extraction and the mixture density output. Additionally we train the model for real time processing using zero padding. In comparison to TEAM-LM we replace the transformer with a maximum pooling operation and remove the event token.
We evaluate two different representations for the position encoding. In the first, we concatenated the positions to the feature vectors as proposed by van den Ende & Ampuero (2020). In the second, we add the position embeddings element-wise to the feature vectors as for TEAM-LM. In both cases, we run a three-layer perceptron over the combined feature and position vector for each station, before applying the pooling operation.
We use the fast magnitude estimation approach (Kuyuk & Allen, 2013) as a classical, i.e., non deep-learning, baseline for magnitude. The magnitude is estimated from the horizontal peak displacement in the first seconds of the P wave. As this approach needs to know the hypocentral distance, it requires knowledge of the event location. We simply provide the method with the catalog hypocenter. While this would not be possible in real time, and therefore gives the method an unfair advantage over the deep learning approaches, it allows us to focus on the magnitude estimation capabilities. Furthermore, in particular for Italy and Japan, the high station density usually allows for sufficiently well constrained location estimates at early times. For a full description of this baseline, see supplement section SM 2.
As a classical location baseline we employ NonLinLoc (Lomax et al., 2000) with the 1D velocity models from Graeber & Asch (1999) (Chile), Ueno et al. (2002) (Japan) and Matrullo et al. (2013) (Italy). For the earliest times after the event detection usually only few picks picks are available. Therefore we apply two heuristics. Until at least 3/5/5 (Chile/Japan/Italy) picks are available, the epicenter is estimated as the arithmetic mean of the stations with picked arrivals so far, while the depth is set to the median depth in the training data set. Until at least 4/7/7 picks are available, we apply NonLinLoc, but fix the depth to the median depth in the data set. We require higher numbers of picks for Italy and Japan, as the pick quality is lower than in Chile but the station density is higher. This leads to worse early NonLinLoc estimates in Italy and Japan compared to Chile, but improves the performance of the heuristics.

Magnitude estimation performance
We first compare the estimation capabilities of TEAM-LM to the baselines in terms of magnitude ( Figure 3). We evaluate the models at fixed times t = 0.5 s, 1 s, 2 s, 4 s, 8 s, 16 s, 25 s after the first P arrival at any station in the network. In addition to presenting selected results here, we provide tables with the results of further experiments in the supplementary material (Tables SM 6-SM 16).
TEAM-LM outperforms the classical magnitude baseline consistently. On two data sets, Chile and Italy, the performance of TEAM-LM with only 0.5 s of data is superior to the baseline with 25 s of data. Even on the third data set, Japan, TEAM-LM requires only approximately a quarter of the time to reach the same precision as the classical baseline and achieves significantly higher precision after 25 s. The RMSE for TEAM-LM stabilizes after 16 s for all data sets Note that the uncertainty of the provided means is by a factor √ 6 smaller than the given standard deviation, due to the number of samples. We provide no standard deviation for the baseline, as it does not depend on a model initialization. Figure 4: RMSE comparison of the TEAM-LM mean magnitude predictions for different magnitude buckets. Linestyles indicate the model type: trained only on the target data (solid line), using transfer learning (dashed), classical baseline (dotted). For Chile/Italy/Japan we count events as small if their magnitude is below 3.5/3.5/4 and as large if their magnitude is at least 5.5/5/6. The time indicates the time since the first P arrival at any station, the RMSE is provided in magnitude units [m.u.].
with final values of 0.08 m.u. for Chile, 0.20 m.u. for Italy and 0.22 m.u. for Japan. The performance differences between TEAM-LM and the classical baseline result from the simplified modelling assumptions for the baseline. While the relationship between early peak displacement and magnitude only holds approximately, TEAM-LM can extract more nuanced features from the waveform. In addition, the relationship for the baseline was originally calibrated for a moment magnitude scale. While all magnitude scales have an approximate 1:1 relationship with moment magnitude, this might introduce further errors.
We further note that the performance of the classical baseline for Italy are consistent with the results reported by Festa et al. (2018). They analyzed early warning performance in a slightly different setting, looking only at the 9 largest events in the 2016 Central Italy sequence. However, they report a RMSE of 0.28 m.u. for the PRESTO system 4 s after the first alert, which matches approximately the 8 s value in our analysis. Similarly, Leyton et al. (2018) analyze how fast magnitudes can be estimated in subductions zones, and obtain values of 0.01 ± 0.28 across all events and −0.70 ± 0.30 for the largest events (M w > 7.5) at 30 s after origin time. This matches the observed performance of the classical baseline for Japan. For Chile, our classical baseline performs considerably worse, likely caused by the many small events with bad SNR compared to the event set considered by Leyton et al. (2018). However, TEAM-LM still outperforms the performance numbers reported by Leyton et al. (2018) by a factor of more than 2.
Figure 5: True and predicted magnitudes without upsampling or transfer learning (left column), with upsampling but without transfer learning (middle column) and with upsampling and transfer learning (right column). All plots show predictions after 8 seconds. In the transfer column for Chile and Japan we show results after fine-tuning on the target data set; for Italy we show results from the model without fine-tuning as this model performed better. For the largest events in Italy (M > 4.5) we additionally show the results after fine-tuning with pale blue dots. We suspect the degraded performance in the fine tuned model results from the fact, that the largest training event (M W = 6.1) is considerably smaller than the largest test event (M W = 6.5). Vertical lines indicate the borders between small, medium and large events as defined in Figure 4. The orange lines show the running 5th, 50th and 95th percentile in 0.2 m.u. buckets. Percentile lines are only shown if sufficiently many data points are available. The very strong outlier for Japan (true ∼7.3, predicted ∼3.3) is an event far offshore (>2000 km).
Improvements for TEAM-LM in comparison to the deep learning baseline variants are much smaller than to the classical approach. Still, for the Japan data set at late times, TEAM-LM offers improvements of up to 27 % for magnitude. For the Italy data set, the baseline variants are on par with TEAM-LM. For Chile, only the baseline with position embeddings is on par with TEAM-LM. Notably, for the Italy and Japan data sets, the standard deviation between multiple runs with different random model initialization is considerably higher for the baselines than for TEAM-LM ( Figure 3, error bars). This indicates that the training of TEAM-LM is more stable with regard to model initialization.
The gains of TEAM-LM can be attributed to two differences: the transformer for station aggregation and the position embeddings. In our experiments we ruled out further differences, e.g. size and structure of the feature extraction CNN, by using identical network architectures for all parts except the feature combination across stations. We think that gains from using a transformer can be explained with its attention mechanism. The attention allows the transformer to focus on specific stations, for example the stations which have recorded the longest waveforms so far. In contrast, the maximum pooling operation is less flexible. We suspect that the high gains for Japan result from the wide spatial distribution of seismicity and therefore very variable station distribution. While in Italy most events are in Central Italy and in Chile the number of stations are limited, the seismicity in Japan occurs along the whole subduction zone with additional onshore events. This complexity can likely be handled better with the flexibility of the transformer than using a pooling operation.
In many use cases, the performance of magnitude estimation algorithms for large magnitude events is of particular importance. In Figure 4 we compare the RMSE of TEAM-LM and the classical baselines for subsets binned by catalog magnitude and observe a clear dependence on the event magnitude. For all data sets the RMSE for large events is higher than for intermediate sized events, which is again higher than for small events. On the other hand the decrease in RMSE over time is strongest for larger events. This general pattern can also be observed for the classical baseline, even though the difference in RMSE between magnitude buckets is smaller. As both variants of the deep learning baseline show very similar trends to TEAM-LM, we omit them from this discussion.
We discuss two possible causes for these effects: (i) the magnitude distribution in the training set restricts the quality of the model optimization, (ii) inherent characteristics of large events. Cause (i) arise from the Gutenberg-Richter distribution of magnitudes. As large magnitudes are rare, the model has significantly less examples to learn from for large magnitudes than for small ones. This should impact the deep learning models the strongest, due to their high number of parameters. Cause (ii) has a geophysical origin. As large events have longer rupture durations, the information gain from longer waveform recordings is larger for large events. We probe the likely individual contributions of these causes in the following.
Estimations for large events not only show lower precision, but are also biased ( Figure 5, middle column). For Chile and Italy a clear saturation sets in for large events. Interestingly the saturation starts at different magnitudes, which are around 5.5 for Italy and 6.0 for Chile. For Japan, events up to magnitude 7 are predicted without obvious bias. This saturation behavior is not only visible for TEAM-LM, but has also been observed in prior studies, e.g., in Mousavi & Beroza (2020) (their Fig. 3, 4). In their work, with a network trained on significantly smaller events, the saturation already set in around magnitude 3. The different saturation thresholds indicate that the primary cause for saturation is not the longer rupture duration of large events or other inherent event properties, as in cause (ii), but is instead likely related to the low number of training examples for large events, rendering it nearly impossible to learn their general characteristics, as in cause (i). This explanation is consistent with the much higher saturation threshold for the Japanese data set, where the training data set contains a comparably large number of large events, encompassing the year 2011 with the Tohoku event and its aftershocks.
As a further check of cause (i), we trained models without upsampling large magnitude events during training, thereby reducing the occurence of large magnitude events to the natural distribution observed in the catalog ( Figure 5, left column). While the overall performance stays similar, the performance for large events is degraded on each of the data sets. Large events are on average underestimated even more strongly. We tried different upsampling rates, but were not able to achieve significantly better performance for large events than the configuration of the preferred model presented in the paper. This shows that upsampling yields improvements, but can not solve the issue completely, as it does not introduce actual additional data. On the other hand, the performance gains for large events from upsampling seem to cause no observable performance drop for smaller event. As the magnitude distribution in most regions approximately follows a Gutenberg-Richter law with b ≈ 1, upsampling rates similar to the ones used in this paper will likely work for other regions as well.

Location estimation performance
We evaluate the epicentral error distributions in terms of the 50 th , 90 th , 95 th and 99 th error percentiles ( Figure 6). In terms of the median epicentral error, TEAM-LM outperforms all baselines in all cases, except for the classical baseline at late times in Italy. For all data sets, TEAM-LM shows a clear decrease in median epicentral error over time. The decrease is strongest for Chile, going from 19 km at 0.5 s to 2 km at 25 s. For Italy the decrease is from 7 km to 2 km, for Japan from 22 km to 14 km. For all data sets the error distributions are heavy tailed. While for Chile even the errors at high quantiles decrease considerably over time, these quantiles stay nearly constant for Italy and Japan.
Similar to the difficulties for large magnitudes, the characteristics of the location estimation point to insufficient training data as source of errors. The Chile data set covers the smallest region and has by far the lowest magnitude of completeness, leading to the highest event density. Consequently the location estimation performance is best and outliers are very rare. For the Italy and Japan data sets, significantly more events occurred in regions with only few training events, causing strong outliers. The errors for the Japanese data set are highest, presumably related to the large number of offshore events with consequently poor azimuthal coverage.
We expect a further difference from the number of unique stations. While for a small number of unique stations, as in the Chile data set, the network can mostly learn to identify the stations using their position embeddings, it might be unable to do so for a larger number of stations with fewer training examples per station. Therefore the task is significantly more complicated for Italy and Japan, where the concept of station locations has to be learned simultaneously to the localization task. This holds true even though we encode the station locations using continuously varying position embeddings. Furthermore, whereas for moderate and large events waveforms from all stations of the Chilean network will show the earthquake and can contribute information, the limitation to 25 stations of the current TEAM-LM implementation does not allow a full exploitation of the information contained in the hundreds of recordings of larger events in the Japanese and Italian data sets. This will matter in particular for out-of-network events, where the wavefront curvature and thus event distance can only be estimated properly by considering stations with later arrivals.
Looking at the classical baseline, we see that it performs considerably worse than TEAM-LM in the Chile data set in all location quantiles, better than TEAM-LM in all but the highest quantiles at late times in the Italy data set, and worse than TEAM-LM at late times in the Japan data set. This strongly different behavior can largely be explained with the pick quality and the station density in the different data sets. While the Chile data set contains high quality automatic picks, obtained using the MPX picker (Aldersons, 2004), the Italy data set uses a simple STA/LTA and the Japan data set uses triggers from KiKNet. This reduces location quality for Italy and Japan, in particular in the case of a low number of picks available for location. On the other hand, the very good median performance of the classical approach for Italy can be explained from the very high station density, giving a strong prior on the location. An epicentral error of around 2 km after 8 s is furthermore consistent with the results from Festa et al. (2018). Considering the reduction in error due to the high station density in Italy, we note that the wide station spacing in Chile likely caused higher location errors than would be achievable with a denser seismic network designed for early warning.
In addition to the pick quality, the assumption of a 1D velocity model for NonLinLoc introduces a systematic error into the localization, in particular for the subduction regions in Japan and Chile where the 3D structure deviates considerably from the 1D model. Because of these limitations the classical baseline could be improved by employing more proficient pickers or fine-tuned velocity models. Nonetheless, in particular the results from Chile, where the classical baseline has access to high quality P-picks, suggest that TEAM-LM can, given sufficient training data, outperform classical real-time localization algorithms.
For magnitude estimation no consistent performance differences between the baseline approach with position embeddings and the approach with concatenated coordinates, as originally proposed by van den Ende & Ampuero (2020), are visible. In contrast, for location estimation, the approach with embeddings consistently outperforms the approach with concatenated coordinates. We speculate that the positional embeddings might show better performance because they explicity encode information on how to interpolate between locations at different scales, enabling an improved exploitation of the information from stations with few or no training examples. This is more important for location estimation, where an explicit notion of relative position is required. In contrast, magnitude estimation can use further information, like frequency content, which is less position dependent.

Transfer learning
A common strategy for mitigating data sparsity is the injection of additional information from related data sets through transfer learning (Pan & Yang, 2009), in our use case waveforms from other source regions. This way the model is supposed to be taught the properties of earthquakes that are consistent across regions, e.g., attenuation due to geometric spreading or the magnitude dependence of source spectra. Note that a similar knowledge transfer implicitly is part of the classical baseline, as it was calibrated using records from multiple regions.
Here, we conduct a transfer learning experiment inspired by the transfer learning used for TEAM. We first train a model jointly on all data sets and then fine-tune it to each of the target data sets. This way, the model has more training examples, which is of special relevance for the rare large events, but still is adapted specifically to the target data set.
As the Japan and Italy data sets contain acceleration traces, while the Chile data set contains velocity traces, we first integrate the Japan and Italy waveforms to obtain velocity traces. This does not have a significant impact on the model performance, as visible in the full results tables (Tables SM 6 to SM 9).
Transfer learning reduces the saturation for large magnitudes ( Figure 5, right column). For Italy the saturation is even completely eliminated. For Chile, while the largest magnitudes are still underestimated, we see a clearly lower level of underestimation than without transfer learning. Results for Japan for the largest events show nearly no difference, which is expected as the Japan data set contains the majority of large events and therefore does not gain significant additional large training examples using transfer learning. The positive impact of transfer learning is also reflected in the lower RMSE for large and intermediate events for Italy and Chile (Figure 4). These results do not only offer a way of mitigating saturation for large events, but also represent further evidence for data sparsity as the reason for the underestimation.
We tried the same transfer learning scheme for mitigating mislocations (Figure 6). For this experiment we shifted the coordinates of stations and events such that the data sets spatially overlap. We note that this shifting is not expected to have any influence on the single data set performance, as the relative locations of events and stations within a data set stay unchanged and nowhere the model uses absolute locations. The transfer learning approach is reasonable, as mislocations might result from data sparsity, similarly to the underestimation of large magnitudes. However, none of the models shows significantly better performance than the preferred models, and in some instances performance even degrades. We conducted additional experiments where shifts were applied separately for each event, but observed even worse performance.
We hypothesize that this behaviour indicates that the TEAM-LM localization does not primarily rely on travel time analysis, but rather employs some form of fingerprinting of earthquakes. These fingerprints could be specific scattering patterns for certain source regions and receivers. Note that similar fingerprints are exploited in the traditional template matching approaches (e.g. Shelly et al. (2007)). While the travel time analysis should be mostly invariant to shifts and therefore be transferable between data sets, the fingerprinting is not invariant to shifts. This would also explain why the transfer learning, where all training samples were already in the pretraining data set and therefore their fingerprints could be extracted, outperforms the shifting of single events, where fingerprints do not relate to earthquake locations.

Multi-task learning
Another common method to improve the quality of machine learning systems in face of data sparsity is multi-task learning (Ruder, 2017), i.e., having a network with multiple outputs for different objectives and training it simultaneously on all objectives. This approach has previously been employed for seismic source characterisation (Lomax et al., 2019), but without an empirical analysis on the specific effects of multi-task learning.
We perform an experiment, in which we train TEAM-LM to predict magnitude and location concurrently. The feature extraction and the transformer parts are shared and only the final MLPs and the mixture density networks are specific to the task. This method is known as hard parameter sharing. The intuition is that the individual tasks share some similarity, e.g., in our case the correct estimation of the magnitude likely requires an assessment of the attenuation and geometric spreading of the waves and therefore some understanding of the source location. This similarity is then expected to drive the model towards learning a solution for the problem that is more general, rather than specific to the training data. The reduced number of free parameters implied by hard parameter sharing is also expected to improve the generality of the derived model, if the remaining degrees of freedom are still sufficient to extract the relevant information from the training data for each sub-task.
Unfortunately, we actually experience a moderate degradation of performance for either location or magnitude in any data set (Tables SM 6 to SM 12) when following a multi-task learning strategy. The RMSE of the mean epicenter estimate increases by at least one third for all times and data sets, and the RMSE for magnitude stays nearly unchanged for the Chile and Japan data sets, but increases by ∼20% for the Italy data set. Our results therefore exhibit a case of negative transfer.
While it is generally not known, under which circumstances multi-task learning shows positive or negative influence (Ruder, 2017), a negative transfer usually seems to be caused by insufficiently related tasks. In our case we suspect that while the tasks are related in a sense of the underlying physics, the training data set is large enough that similarities relevant for both tasks can be learned already from a single objective. At the same time, the particularities of the two objectives can be learned less well. Furthermore, we earlier discussed that both magnitude and location might not actually use travel time or attenuation based approaches, but rather frequency characteristics for magnitude and a fingerprinting scheme for location. These approaches would be less transferable between the two tasks. We conclude that hard parameter sharing does not improve magnitude and location estimation. Future work is required to see if other multi-task learning schemes can be applied beneficially.

Location outlier analysis
As all location error distributions are heavy tailed, we visually inspect the largest deviations between predicted and catalog locations to understand the behavior of the localization mechanism of TEAM-LM. We base this analysis on the Chile data set (Figure 7), as it has generally the best location estimation performance, but observations are similar for the other data sets ( Figures SM 5 and SM 6).
Nearly all mislocated events are outside the seismic network and location predictions are generally biased towards the network. This matches the expected errors for traditional localization algorithms. In contrast to traditional algorithms, events are not only predicted to be closer to the network, but they are also predicted as lying in regions with a higher event density in the training set (Figure 7, inset). This suggests that not enough similar events were included in the training data set. Similarly, Kriegerowski et al. (2019) observed a clustering tendency when predicting the location of swarm earthquakes with deep learning.
We investigated two subgroups of mislocated events: the Iquique sequence, consisting of the Iquique mainshock, foreshocks and aftershocks, and mine blasts. The Iquique sequence is visible in the north-western part of the study area. All events are predicted approximately 0.5 • too far east. The area is both outside the seismic network and has no events in the training set. This systematic mislocation may pose a serious threat in applications, such as early warning, when confronted with a major change in the seismicity pattern, as is common in the wake of major earthquakes or during sudden swarm activity, which are also periods of heightened seismic hazard.
For mine blasts, we see one mine in the northeast and one in the southwest (marked by red circles in Figure 7). While all events are located close by, the location are both systematically mispredicted in the direction of the network and exhibit scatter. Mine-blasts show a generally lower location quality in the test set. While they make up only ∼1.8% of the test set, they make up 8% of the top 500 mislocated events. This is surprising as they occur not only in the test set, but also in similar quantities in the training set. We therefore suspect that the difficulties are caused by the strongly different waveforms of mine blasts compared to earthquakes. One waveform of each a mine blast and an earthquake, recorded at similar distances are shown as inset in Figure 7. While for the earthquake both a P and S wave are visible, the S wave can not be identified for the mine blast. In addition, the mine blast exhibits a strong surface wave, which is not visible for the earthquake. The algorithm therefore can not use the same features as for earthquakes to constrain the distance to a mine blast event.

The impact of data set size and composition
Our analysis so far showed the importance of the amount of training data. To quantify the impact of data availability on magnitude and location estimation, we trained models only using fractions of the training and validation data ( Figure  8). We use the Chile data set for this analysis, as it contains by far the most events. We subsample the events by only using each k th event in chronological order, with k = 2, 4, 8, 16, 32, 64. This strategy approximately maintains the magnitude and location distribution of the full set. We point out, that TEAM-LM only uses information of the event under consideration and does not take the events before or afterwards into account. Therefore, the 'gaps' between events introduced by the subsampling do not negatively influence TEAM-LM.
For all times after the first P arrival, we see a clear increase in the magnitude-RMSE for a reduction in the number of training samples. While the impact of reducing the data set by half is relatively small, using only a quarter of the data already leads to a twofold increase in RMSE at late times. Even more relevant in an early warning context, a fourfold smaller data sets results in an approximately fourfold increase in the time needed to reach the same precision as with the full data. This relationship seems to hold approximately across all subsampled data sets: reducing the data set k fold increases the time to reach a certain precision by a factor of k.
We make three further observations from comparing the predictions to the true values ( Figure SM 7). First, for nearly all models the RMSE changes only marginally between 16 s and 25 s, but the RMSE of this plateau increases significantly with a decreasing number of training events. Second, the lower the amount of training data, the lower is the saturation threshold above which all events are strongly underestimated. In addition, for 1/32 and 1/64 of the full data set, an 'inverse saturation' effect is noticeable for the smallest magnitudes. Third, while for the full data set and the largest subsets all large events are estimated at approximately the saturation threshold, if at most one quarter of the training data is used, the largest events even fall significantly below the saturation threshold. For the models trained on the smallest subsets (1/8 to 1/64), the higher the true magnitude the lower the predicted magnitude becomes. We assume that the larger the event is, the further away from the training distribution it is and therefore it is estimated approximately at the  All models were evaluated on the full Chilean test set. We note that the variance of the curves with fewer data is higher, due to the increased stochasticity from model training and initialization. most dense region of the training label distribution. These observations support the hypothesis that underestimations of large magnitudes for the full data set are caused primarily by insufficient training data.
While the RMSE for epicenter estimation shows a similar behavior as the RMSE for magnitude, there are subtle differences. If the amount of training data is halved, the performance only degrades mildly and only at later times. However, the performance degradation is much more severe than for magnitude if only a quarter or less of the training data are available. This demonstrates that location estimation with high accuracy requires catalogs with a high event density.
The strong degradation further suggests insights into the inner working of TEAM-LM. Classically, localization should be a task where interpolation leads to good results, i.e., the travel times for an event in the middle of two others should be approximately the average between the travel times for the other events. Following this argument, if the network would be able to use interpolation, it should not suffer such significant degradation when faced with fewer data. This provides further evidence that the network does not actually learn some form of triangulation, but only an elaborate fingerprinting scheme, backing the finding from the qualitative analysis of location errors.

Training TEAM-LM on large events only
Often, large events are of the greatest concerns, and as discussed, generally showed poorer performance because they are not well represented in the training data. It therefore appears plausible that a model optimized for large events might perform better than a model trained on both large and small events. In order to test this hypothesis, we employed an extreme version of the upscaling strategy by training a set of models only on large events, which might avoid tuning the model to seemingly irrelevant small events. In fact, these models perform significantly worse than the models trained on the full data set, even for the large events (Tables SM 6 to SM 12). Therefore even if the events of interest are only the large ones, training on more complete catalogs is still beneficial, presumably by giving the network more comprehensive information on the regional propagation characteristics and possibly site effects.

Interpretation of predicted uncertainties
So far we only analyzed the mean predictions of TEAM-LM. As for many application scenarios, for example early warning, quantified uncertainties are required, TEAM-LM outputs not only these mean predictions, but a probability density. Figure 9 shows the development of magnitude uncertainties for events from different magnitude classes in the Chile data set. As we average over multiple events, each set of lines can be seen as a prototype event of a certain magnitude. Figure 9: Magnitude predictions and uncertainties in the Chile data set as a function of time since the first P arrival. Solid lines indicate median predictions, while dashed lines show 20th and 80th quantiles of the prediction. Each color represents a different magnitude bucket. For each magnitude bucket, we sampled 1,000 events around this magnitude and combined their predictions. If less than 1,000 events were available within ±0.5 m.u. of the bucket center, we use all events within this range. We only use events from the test set. To ensure that the actual uncertainty distribution is visualized, rather than the distribution of magnitudes around the bucket center, each prediction is shifted by the magnitude difference between bucket center and catalog magnitude.
For all magnitude classes the estimation shows a sharp jump at t = 0, followed by a slow convergence to the final magnitude estimate. We suspect that the magnitude estimation always converges from below, as due to the Gutenberg-Richter distribution, lower magnitudes are more likely a priori. The uncertainties are largest directly after t = 0 and subsequently decrease, with the highest uncertainties for the largest events. As we do not use transfer learning in this approach, there is a consistent underestimation of the largest magnitude events, visible from the incorrect median predictions for magnitudes 5 and 6.
While the Gaussian mixture model is designed to output uncertainties, it cannot be assumed that the predicted uncertainties are indeed well calibrated, i.e., that they actually match the real error distribution. Having well calibrated uncertainties is crucial for downstream tasks that rely on the uncertainties. Neural networks trained with a log-likelihood loss generally tend to be overconfident (Snoek et al., 2019;Guo et al., 2017), i.e., underestimate the uncertainties. This overconfidence is probably caused by the strong overparametrization of neural network models. To assess the quality of our uncertainty estimations for magnitude, we use the observation that for a specific event i, the predicted Gaussian mixture implies a cumulative distribution function F i pred : R → [0, 1]. Given the observed magnitude y i true , we can calculate u i = F i pred (y i true ). If y i true is indeed distributed according to F i pred , then u i needs to be uniformly distributed on [0, 1]. We test this based on the u i of all events in the test set using P-P plots ( Figure 10). Further details on the method can be found in the supplementary material (Section SM 3). Note that good calibration is a necessary but not sufficient condition for a good probabilistic forecast. An example of a perfectly calibrated but mostly useless probabilistic prediction would be the marginal probability of the labels. Figure 10 shows the P-P plots of u in comparison to a uniform distribution. For all data sets and all times the model is signficantly miscalibrated, as estimated using Kolmgorov-Smirnov test statistics (Section SM 3). Miscalibration is considerably stronger for Italy and Japan than for Chile. More precisely, the model is always overconfident, i.e., estimates narrower confidence bands than the actually observed errors. Further, in particular at later times, the model is biased towards underestimating the magnitudes. This is least visible for Chile. We speculate that this is a result of the large training data set for Chile, which ensures that for most events the density of training events in their magnitude range is high.
To mitigate the miscalibration, we trained ensembles (Hansen & Salamon, 1990), a classical method to improve calibration. Instead of training a single neural network, a set of n neural networks, in our case n = 10, are trained, which all have the same structure, but different initialization and batching in training. The networks therefore represent a sample of size n from the posterior distribution of the model parameters given the training data. For Italy and Japan, this improves calibration considerably (Figure 10). For Chile, the ensemble model, in contrast to the single model, exhibits underconfidence, i.e., estimates too broad uncertainty bands.
The maximum distance between the empirical cumulative distribution function of u and a uniformly distributed variable d ∞ is the test statistic of the Kolmogorov-Smirnov test. While d ∞ is reduced by nearly half for some of the ensemble Figure 10: P-P plots of the CDFs of the empirical quantile compared to the expected uniform distribution. The P-P plot shows (CDF ui (z), CDF uniform (z)) for z ∈ [0, 1]. The expected uniform distribution is shown as the diagonal line, the misfit is indicated as shaded area. The value in the upper corner provides d ∞ , the maximum distance between the diagonal and the observed CDF. d ∞ can be interpreted as the test statistic for a Kolmogorov-Smirnov test. Curves consistently above the diagonal indicate a bias to underestimation, and below the diagonal to overestimation. Sigmoidal curves indicate over-confidence, mirrored sigmoids indicate under-confidence. See supplementary section SM 3 for a further discussion of the plotting methodology and its connection to the Kolmogorov-Smirnov test. Figure 11: The figure shows 90th percent confidence areas for sample events around 5 example locations. For each location the 5 closest events are shown. Confidence areas belonging to the same location are visualized using the same color. Confidence areas were chosen as curves of constant likelihood, such that the probability mass above the likelihood equals 0.9. To visualize the result in 2D we marginalize out the depth. Triangles denote station locations for orientation. The top row plots show results from a single model, while the bottom row plots show results from an ensemble of 10 models. results, the Kolmogorov-Smirnov test indicates, that even the distributions from the ensemble models deviate highly significantly from a uniform distribution (p 10 −5 ). A table with d ∞ for all experiments can be found in the supplementary material (Table SM 9).
To evaluate the location uncertainties qualitatively, we plot confidence ellipses for a set of events in Chile (Figure 11). Again we compare the predictions from a single model to the predictions of an ensemble. At early times, the uncertainty regions mirror the seismicity around the station with the first arrival, showing that the model correctly learned the prior distribution. Uncertainty ellipses at late times approximately match the expected uncertainty ellipses for classical methods, i.e., they are small and fairly round for events inside the seismic network, where there is good azimuthal coverage, and larger and elliptical for events outside the network. Location uncertainties are not symmetric around the mean prediction, but show higher likelihood towards the network than further outwards. Location errors for the ensemble model are more smooth than from the single model, but show the same features. The uncertainty ellipses are slightly larger, suggesting that the single model is again overconfident.
In addition to improving calibration, ensembles also lead to slight improvements regarding the accuracy of the mean predictions (Tables SM 6 to SM 12). Improvements in terms of magnitude RMSE range up to 10%, for epicentral location error up to 20%. Due to the high computational demand of training ensembles, all other results reported in this paper are calculated without ensembling.

Conclusion
In this study we adapted TEAM to build TEAM-LM, a real time earthquake source characterization model, and used it to study the pitfalls and particularities of deep learning for this task. We showed that TEAM-LM achieves state of the art in magnitude estimation, outperforming both a classical baseline and a deep learning baseline. Given sufficiently large catalogs, magnitude can be assessed with a standard deviation of ∼0.2 magnitude units within 2 s of the first P arrival and a standard deviation of 0.07 m.u. within the first 25 s. For location estimation, TEAM-LM outperforms a state of the art deep learning baseline and compares favorably with a classical baseline.
Our analysis showed that the quality of model predictions depends crucially on the training data. While performance in regions with abundant data is excellent, in regions of data sparsity, prediction quality degrades significantly. For magnitude estimation this effect results in the underestimation of large magnitude events; for location estimation events in regions with few or no training events tend to be mislocated most severely. This results in a heavy tailed error distribution for location estimation. Large deviations in both magnitude and location estimation can have significant impact in application scenarios, e.g., for early warning where large magnitudes are of the biggest interest.
Following our analysis, we propose a set of best practices for building models for fast earthquake source characterisation: 1. Build a comprehensive evaluation platform. Put a special focus on outliers and rare or large events. Analyze which impact outliers or out of distribution events will have for the proposed application.
2. Use very large training catalogs, spanning long time spans and having a low magnitude of completeness. If possible, employ transfer learning. We hope the catalogs used in this study can give a starting point for transfer learning.
3. Use training data augmentation, especially upsampling of large events, which improves prediction performance in face of label sparsity at virtually no cost.
4. If probabilistic estimates are required, use deep ensembles to improve the model calibration.
5. When using deep learning for location estimation, put special emphasis on monitoring possible distribution shifts between training data and application.
While these points give guidance for training current models they also point to further directions for methodological advances. First, our transfer learning scheme is fairly simple. More refined and targeted schemes could increase the amount of information sharable across data sets considerably. We further expect major improvements from training with simulated data, but are aware that generating realistic, synthetic seismogramms, especially for large events, poses major challenges. Another promising alternative might be to move away from the paradigm of black box modeling, i.e., training algorithms that are built solely by fitting recorded data. Instead, incorporation of physical knowledge and a move towards physics informed deep learning methods seems promising (Raissi et al., 2019). However, physics informed neural networks are still in their infancy and the application to seismic tasks still needs to be developed.

SM 1 Data sources SM 2 Classical magnitude estimation baseline
For magnitude estimation we compare TEAM-LM to a classical baseline. To this end we use the peak displacement based approach proposed by Kuyuk & Allen (2013). At each station, we bandpass filter the signal between 0.5 Hz and 3 Hz and discard traces with insufficient signal to noise ratio. We extract peak displacement P D from the horizontal components in the first 6 s of the P wave, while only including samples before the S onset. We use the relationship from Kuyuk & Allen (2013) to estimate magnitudes from peak displacement. We use c 1 = 1.23, c 2 = 1.38 and σ = 0.31 from Kuyuk & Allen (2013). These parameters were calibrated using data from California and Japan, but the authors state that the relationship can be applied to earthquake source zones around the world. To account for a constant offset between different magnitude scales, we optimized c 3 separately for each dataset such that the predictions do not have a systematic bias compared to the ground truth.
We average the predictions from multiple stations, effectively assuming independence between the predictions. To obtain earlier predictions, we already calculate magnitude estimates at a station once at least 1 s of P wave data has been recorded. We assign higher weights to stations with longer P wave records, with weights linearly increasing from 0.11 for 1 s of waveforms, to 1.0 for 6 s of data. Thereby, while getting early estimates from the first stations, new data from later stations does not perturb the prediction strongly until enough data has been recorded.
As the estimation relies on the hypocentral distance R between station and event, the method requires an estimate of the hypocentral location. We provide the method with the cataloged hypocentral location. While this is an unrealistically optimistic assumption for an actual real time determination, it allows us to put our focus on the magnitude estimation capabilities. We note that this gives the baseline an advantage compared to TEAM-LM, which has no information on the earthquake location.
For some events in the Chile catalog, the SNR criterion is not fulfilled at any station due to the inclusion of smaller magnitude events and the higher distances between stations and events. For these events the baseline does not issue a magnitude estimation. We exclude these events from the evaluation of the baseline, leading to an optimistic assessment of the performance of the baseline.

SM 3 Calibration estimation
Calibration of a model describes whether the predicted uncertainties match the observed values, i.e., if the observation y true was drawn from a distribution with cumulative distribution function (CDF) F pred . Unfortunately, for each event i, only one prediction observation of the magnitude y i true and one prediction F i pred is available. To this end, we define the random variable u i = F i pred (y i true ). If y i true is distributed according to F i pred than u i must be uniformly distributed on [0,1]. This follows from the definition of the CDF. If F is a CDF and U a uniform random variable on [0, 1], than F −1 (U ) is distributed according to F . We take the u i of all events as samples of a random variable U and compare U to a uniform random variable on [0, 1]. The maximum difference between the empirical CDF of U and a uniform variable is the test statistic of a Kolmogorov-Smirnov test, d ∞ . As the number of events n is large, critical values d α to a confidence threshold α can be estimated as: For α = 10 −5 , this gives values d α of 0.015 (Chile), 0.054 (Japan) and 0.039 (Italy). This is considerably below the observed values d ∞ , even using ensembles, indicating that U differs highly significantly from a uniform distribution.     The event density is calculated using a Gaussian kernel density estimation and does not take into account the event depth. The inset shows the event density at the true event location in comparison to the event density at the predicted event location. The event density is calculated using a Gaussian kernel density estimation and does not take into account the event depth. The inset shows the event density at the true event location in comparison to the event density at the predicted event location. All models were evaluated on the full test dataset.    Table SM 13: Test set root squared mean for hypocentral error. We note that only 4 out of 10 models for the Italy location ensemble converged. We used only the converged models for the ensemble evaluation.