Citation: Madi, M., T. Hara (2023), Determination of hypocenters and magnitudes of local earthquakes
around Comoros, Synopsis of IISEE-GRIPS Master’s Thesis.

DETERMINATION OF HYPOCENTERS ANDMAGNITUDES OF LOCAL EARTHQUAKES AROUNDCOMOROSMariama MADI1
MEE22703

Supervisor: Tatsuhiko HARA2*

ABSTRACT
We have conducted hypocenter determinations and calculated local magnitudes for five events, which
occurred on November 23, 24, and 29, 2019, May 29, 2020, and October 13, 2022, around Comoro
Islands. We retrieved data recorded by the local seismic network and picked P and S arrival times using
the seismograms and the probability density functions of their amplitudes. We determined hypocenters
using two velocity models, which are a global model and a local model, respectively. The focal depth
was set to 10 km considering the difficulty in its determination. The RMS (root-mean-square) residuals
for the local model are found to be smaller, which implies that the local model is more appropriate.
Then we calculated local magnitudes by amplitude measurements of synthetic Wood-Anderson
seismograms. Since duration magnitude is the only one calculated at the Karthala Volcano Observatoryof Comoros, the local magnitude will be the second magnitude to be registered in the observatory
bulletin. Although the calculated local magnitudes are comparable to the duration magnitudes, theircorrelation is unclear. Further analyses are necessary to establish their relation.
The methods and analyses used in this study are likely to improve earthquake monitoring at the Karthala
Volcano Observatory of the Comoros. Further analyses using a larger dataset will be useful formonitoring the seismic and volcanic activities in the study region.
Keywords: Hypocenter, Local magnitude, Comoro Islands
1. INTRODUCTION
The Comoro Islands is a chain of four islands located in the entry of the Mozambique Channel between
the southeastern part of Africa and the northwestern side of Madagascar. The magmatism andmorphology of the islands are comparable to those of the volcanism of the Hawaiian Islands which are
known to have various of stages and erosion development processes. The Comoros hotspot is one of thehotspot volcanoes: (e.g., Steinberger, 2000). The Comoro Islands have significant records of basalts andassociated volcanic and plutonic rocks that compose them (e.g., Thompson & Flower, 1971; Flower,
1971, 1972, 1973a,b; Ludden, 1977), while the islands become younger toward the Western side.
The Gde-Comore-Island is younger than the other three archipelago islands, Moheli, Anjouan,
and Mayotte. The active Karthala volcano induces seismic activities recorded in the Gde-Comore
Island. The volcanic activities in 2004, 2005 (Thivet et al., 2022), 2006 (Morin, et al., 2009), 2007(https://doi.org/10.5479/si.GVP.BGVN200701-233010), and the recent 2022 volcano-seismic activity
1

National Data Center (NDC) and Karthala Volcano Observatory (OVK), under the National Center ofDocumentation and Scientific Research (CNDRS), Moroni, Comoros.
2
Professor, International Institute of Seismology and Earthquake Engineering, Building Research Institute.
* Chief examiner1require more attention in terms of volcano monitoring. The improvement of the volcano earthquake
observation capacity is crucial for the Gde-Comore-Island. In this work, we analyzed the local network
data applying the new method using the PDF (Probability Density Function) of seismic wave amplitudes
Oshima and Takenaka (2020, 2022) to pick P and S arrivals. Using HYPOCENTER software program(Lienert et al., 1986; Lienert & Havskov, 1995) installed in SEISAN Software (Havskov & Ottemoller,
1999; Havskov et al., 2020), we determined the hypocenters of five local events, that is, three events in2019 (May 23, 2019, May 24, 2019, and May 29, 2019), one in 2020 (November 24) and one in 2022
(October 13) around the island and with six stations from the local network. We calculate the local
magnitudes using the SEISAN software and compare them to the duration magnitudes obtained fromthe volcano observatory bulletin, which are determined routinely.
2. METHODOLOGYThis chapter explains the methods of picking the P and S arrivals, hypocenter determination, and local
magnitude calculation.
2.1. Phase picking using the probability density function of seismic waveform amplitudes
Oshima and Takenaka (2020) developed a method for picking the P and S arrivals based on theprobability density function (PDF) of the seismic wave amplitudes. Oshima and Takenaka (2022)
developed a new method for enhancing direct waves using the PDF of waveform amplitude. The KLD
Kullback and Leibler, 1951) is used to quantify the dissimilarity between the PDF and the Gaussian orRayleigh distribution.
∞𝑝(𝑥)𝐾𝐿𝐷 (𝑝, 𝑞) = ∫−∞ 𝑝(𝑥)log
𝑑𝑥 ,
(1)
𝑞(𝑥)

where, 𝑝(𝑥) and 𝑞(𝑥) are the PDFs, and the similarity between 𝑝(𝑥) and 𝑞(𝑥) is quantified.
Oshima and Takenaka (2022) calculated the KLD between the PDF of the absolute amplitude and its
corresponding Rayleigh distribution as follows:𝑝(𝑥𝑖 )𝐾𝐿𝐷 (𝑝, 𝑞) = ∑𝑛𝑏𝑖𝑛(2)
𝑖=1 𝑝(𝑥𝑖 )log 𝑞(𝑥 ),𝑖where, nbin is the number of bin, which is √𝑛𝑤𝑖𝑛𝑑; nwind is the number of the normalized absolute
amplitude data, 𝑥(𝑘) (𝑘 = 0, 1, . . . , 𝑛𝑤𝑖𝑛𝑑 − 1); 𝑥𝑖 is the central value of the i-th bin; 𝑝(𝑥𝑖 ) is thenumber of the data in the i-th bin devided by nwind; and 𝑞(𝑥𝑖 ) is the value of the Rayleigh distribution.
The baseline is corrected before calculations. The Rayleigh distribution is represented by𝑥𝑥2(3)
𝑅(𝑥) = 2 exp(− 2),and the parameter σ is calculated as𝜎

2𝜎𝜎 = √∑𝑛𝑤𝑖𝑛𝑑
𝑘=1 {𝑥(𝑘)}²/2𝑛𝑤𝑖𝑛𝑑,(4)

Oshima and Takenaka (2020) used the two KLD indexes. The first index is KLD1 defined as the KLD
between the PDF of the observed amplitude waveform, and the Raleigh or Gaussian distribution. In thisstudy, we used the Rayleigh distribution following Oshima and Takenaka (2022). The second index is
KLD2 defined as the KLD between the PDFs for the adjacent time windows on a single seismic trace.
They used KLD1 for the P-phase picking on the amplitude of the vertical component. Both KLD1 andKLD2 are used for S-phase picking. Although Oshima and Takenaka (2020) showed that this technique
was applicable to automatic phase picking, we used these KLDs to support manual picking.
2.2. Hypocenter Determination
We used the HYPOCENTER software program (Lienert et al., 1986; Lienert & Havskov, 1995) todetermine the hypocenters using the P and S arrival times from the local data. It is included in the
SEISAN software (Havskov & Ottemoller, 1999; Havskov et al., 2020), which can be used for2earthquake data analyses at seismic observations. HYPOCENTER is a default program for hypocenter
determination in the SEISAN, for which a file including information such as station locations andvelocity models should be prepared before analyses.
2.2.1. Velocity Models
We used the two velocity models to determine the hypocenters using HYPOCENTER: the iasp91 model(Kennett and, Engdahl, 1991) and a model set at -11.5°Latitude and 43.5° Longitude in CRUST 1.0(Laske, et al., 2013). We used “Crust 1.0 and Google Earth” (https://github.com/jrleeman/Crust1.0) to
extract the model parameters of CRUST1.0. For the latter, we applied the model parameters in the crustand upper mantle and refer to the model as CR1. Table 1 shows CR1 model. The thickness of ias91modelis 35 km, while for CR1 it is thinner (14.15 km). Model CR1 includes a thin relatively slow velocity
layer at the top.
Table 1. The Model parameters of the model from CRUST 1.0 (Laske, et al., 2013), which we referto as the CR1 model.
Depth (km) to the top of the layer

Vp (km)Vs (km)

0

5

2.7

1.41

6.5

3.7

4.24

7.1

4.05

14.15

8.13

4.51

2.3. Local Magnitude Calculation
We used SEISAN to measure the amplitudes and calculate the local magnitudes. The Synthetic WoodAnderson seismograms are calculated by deconvolution using the instrumental responses and theconvolution of the Wood-Anderson seismograph response. Following Havskov and Ottemoller (2010),
we used vertical component seismograms and applied a high-band pass filter in cases where the seismic
signals are difficult to see on the synthetic Wood-Anderson seismograms. In SEISAN, the effect of thehigh-band pass filter is corrected for amplitude measurements. The local magnitude (ML) is calculatedby the following formular (Borman et al 2013; Borman and Dewey, 2014):
𝑀𝐿 = log10 (𝐴) + 1.11 log10 𝑅 + 0.00189 ∗ 𝑅 − 2.9,

(5)where, A is the maximum amplitude (nm) of the synthetic Wood-Anderson seismogram, and R is the
hypocentral distance (km). The coefficients for the first and second terms are from Hutton and Boore(1987), who obtained their local magnitude formula for Southern California.
3. DATA AND PROCESSING PROCEDURE
3.1. Data and Seismic Network
We analyzed five earthquakes that occurred around Comoros Islands. We selected them by considering
the data availability and their relatively large sizes compared to events in the main island. As mentioned
in the Introduction section, the Karthala Volcano Observatory of Comoros routinely determines theduration magnitudes using Seiscomp3 (Helmholtz Centre Potsdam GFZ German Research Centre forGeosciences and gempa GmbH, 2008). We retrieved the waveform data recorded at six stations (CAB,
DEMB, GRIL, MOIN, NBC, NTC) of the Karthala Volcano Observatory of Comoros.
3.2. Procedure of Data Processing
3.2.1. Picking of P and S Arrival Times3

Figure 1 shows examples of the SAC data and KLD1 and KLD2 results. The event occurred on
November 24, 2020, and the station is DEMB. We used 2s time window for the KLD1 calculation and2s and 4s as time windows for the KLD2 calculation. The first to third traces from the top show theoriginal three-component waveforms; the fourth to sixth traces from the top show their KLD1
waveforms; The seventh to nineth traces from the top denotes the KLD2 waveforms for with a is 2 stime window; and the 10th to 12th traces from the top display the KLD2 waveforms with a 4s time
window. The arrivals of P and S waves can be seen on the KLD seismograms. In some cases, P and S
arrivals are not clear on KLD seismograms.
We calculated KLDs by the procedures explained in
2.1. As Oshima and Takenaka (2020) showed if thecoda part is used for the Rayleigh distribution thechanges of KLD1 will be larger. Also it is desirable
to further examine the effects of time widows andnoise levels. The application of band-pass filtering
will be another subject to be studied in future.

Figure 1. Example of the SAC waveform data, and their KLD1 and KLD2. The event is the one that
occurred on November 24, 2020, and the station is DEMB.
3.2.2. Hypocenter Determination
After picking the P and S arrival times, we put the arrival time data in the SEISAN database andconducted hypocenter determination using HYPOCENTER. The station locations are put in
STATIONX.HYP, which is used by HYPOCENTER. We prepared this file for model iasp91 and CR1,
respectively. The waveform data, which were converted to miniSEED format, were registered in the
SEISAN database.We tried to determine the origin time and the hypocenter, however found that it wasdifficult to determine focal depths due to the limited number of data. Therefore, we decided to fix the
focal depth for each event and set it to 10 km in this study.
3.2.3. Amplitude measurement
After conducting the hypocenter determinations, we proceeded with the amplitude measurements using
the vertical component seismograms. We used the SEISAN software (Havskov & Ottemoller, 1999;
Havskov et al., 2020) to calculate the synthetic Wood-Anderson seismograms and used them for manual
measurements after applying high band pass filter (15-25HZ). As stated in section 2.3, the effect of thehigh-band pass filter is corrected in SEISAN.
4. RESULTS AND DISCUSSION
4.1. Hypocenter Determination
We determined hypocenters for the five events presented in Table 2 using the method and procedure
explained in the previous chapters. Figures 2(a) and (b) show the results of these calculations. The errorsof the origin times and epicenters for model CR1 are smaller than or comparable to those for model
iasp91 (Kennett and Engdahl, 1991). The RMSs for model CR1 are smaller than those for model iasp91.
These imply that CR1 is better than iasp91 for hypocenter determination for events in the study region.
Since the number of the analyzed events are small, further investigations are necessary to find a suitablemodel. The focal depths of some hypocenters in the bulletin of the Volcano Observatory are very deep,
say around 50 km or larger than 100 km. Considering the tectonic environment of the study area, suchdeep events are not likely, and probably they are due to the difficulty to determine focal depths.
Therefore for such a case, it is better and reasonable to set a shallow focal depth.

4Table 2. Hypocenters determined using the CR1 and iasp91 (Kennett and Engdahl, 1991) models forthe five events.
Event

1

2

3

4

5

Velocitymodel
crust1.0

Date

Origin time

LatitudeLat.error

Longitude

Lon.error

RMS

22.23.57.5

Time
error
0.530

2019 May 23

iasp912019 May 23-12.096

14.7

42.3946.6

0.23

22.23.55.4

0.660-12.22212.1

42.562

7.5

crust1.0

0.48

2019 May 24

02.17.46.0

0.600-12.68831.4

41.64716.9

0.3

iasp912019 May 24

02.17.45.1

0.610-12.57926

41.93915

0.48

crust1.0

2019 May 29

17.49.48.4

0.330-12.111

5.6

43.9113.8

0.11

iasp912019 May 29

17.49.47.8

1.17-12.0413.3

43.85910.2

0.58

crust1.0

2020 Nov 24

23.45.16.0

0.290-12.485343.18880.07

iasp912020 Nov 24

23.45.13.6

1.73-12.3315.5

42.94722.2

1.15

crust1.02022 Oct 13

10.11.24.7

0.850-12.27414.8

44.24111.3

0.4

iasp912022 Oct 13

10.11.23.3

1.21-12.0514

44.188

9.4

0.62

(a)

(b)

Figure 2. (a) the epicenters determined using the CR1 model. (b) the epicenters determined using
iaps91 model (Kennett & Engdahl, 1991). The focal depths are fixed to 10km.
4.2. Local magnitude calculation
We calculated ML for the five events using the hypocenters obtained for CR1 model. Figure 3 shows
the comparison between ML and Md calculated by the Karthala Volcano Observatory of the GdeComore Island. The estimates themselves for ML are comparable to those of Md. Also, the estimates of
Md are larger than those of ML for all the events. Considering the difficulty associated with amplitude
measurements for small events, it is necessary to calculate Md and register it into the earthquake bulletin
continuously. Further local magnitude analyses of more events are necessary to evaluate the accuracyof Md and the relation between Md and Ml.

Figure 3. Comparison between ML calculated in this study and Md from the observatory bulletin.
5. CONCLUSIONSIn this study, we determined hypocenters and calculated local magnitudes for five earthquakes, that
occurred in 2019, 2020, and 2022, using waveform data recorded at the stations of the Karthala Volcano

5Observatory. First, we picked P and S arrival times using the recorded seismograms and their KLD
seismograms, which were calculated following the methods of Oshima and Takenaka (2020, 2022).
Then we determined hypocenters using HYPOCENTER program in SEISAN software, in the databaseof which P and S arrival times waveforms are put. We used two velocity models for our hypocenter
determination, which are CR1 from CRUST 1.0 (Laske et al., 2013) and iasp91 (Kennett and Engdahl,1991). We fixed the focal depths at 10km for all the five events. The small RMS residuals imply modelCR1 is better for hypocenter determination.
We measured amplitudes of synthetic Wood-Anderson seismograms calculated for the observed
seismograms using SEISAN. We then calculated ML using the measured amplitudes and the hypocenter
distances from the hypocenters obtained for model CR1; we used the formula of the IASPEI standard
procedure with the coefficients for southern California based on the work of Hutton & Boore (1987).
The calculated ML estimates are comparable to those of Md calculated in the Observatory bulletin. It
seems that the correlation between them is not clear. However, the number of events needs to be
increased to evaluate their relation more precisely.
ACKNOWLEDGEMENTS
I would like to take this opportunity to express my sincere gratitude to Dr. Tatsuhiko Hara forsupervising me during my study and providing me with adequate guidance.
REFERENCES
Bormann, P., & Dewey, J. W. (2014). IS3. 3. NMSOP2. German Research Centre for Geosciences.
Bormann, P., Wendt, S., & DiGiacomo, D. (2013). In New manual of seismological observatory practice
2 (NMSOP2) (pp. 1-259). Deutsches GeoForschungsZentrum GFZ.
Flower, M. (1971). Contributions to Mineralogy and Petrology, 31(4), 335-346.
Flower, M. (1972). Bulletin Volcanologique, 36, 238-250.
Flower, M. (1973a). Chemical Geology, 12(2), 81-98.
Flower, M. (1973b). Contributions to Mineralogy and Petrology, 38(3), 237-260.
Global Volcanism Program. (2007). Bulletin of the Global Volcanism Network, 32:1. Smithsonian
Institution. https://doi.org/10.5479/si.GVP.BGVN200701-233010
Havskov, J., Voss, P. H., & Ottemöller, L. (2020). Seismological Research Letters, 91(3), 1846-1852.
Havskov, J., & Ottemoller, L. (1999). Seismol. Res. Lett, 70(5), 532-534.
Helmholtz Centre Potsdam GFZ German Research Centre for Geosciences and gempa GmbH (2008).
GFZ Data Services. doi:10.5880/GFZ.2.4.2020.003.
Hutton, L. K., & Boore, D. M. (1987). Bulletin of the SSA, 77(6), 2074-2094.
Kennett, B., & Engdahl, E. (1991). Geophysical Journal International, 105(2), 429-465.
Kullback, S., & Leibler, R. A. (1951).The annals of mathematical statistics, 22(1), 79-86.
Laske, G., Masters., G., Ma, Z. & Pasyanos, M. (2013). Abstracts, 15, Abstract EGU2013-2658.
Lienert, B. R., Berg, E., & Frazer, L. N. (1986). HYPOCENTER: Bulletin of the Seismological Societyof America, 76(3), 771-783.
Lienert, B. R., & Havskov, J. (1995). Seismological Research Letters, 66(5), 26-36.
Ludden, J. (1977). Contributions to Mineralogy and Petrology, 64(1), 91-107.
Morin, J., Lavigne, F., Bachelery, P., Finizola, A. and Villeneuve, N., 2009. Shima: The International
Journal of Research into Island Cultures., 3(1), pp.33-53.
Oshima, M., & Takenaka, H. (2020). Bulletin of the Seismological Society of America, 110(2), 763-782.
Oshima, M., & Takenaka, H. (2022). Geophysical Journal International, 231(1), 327-354.
Steinberger, B. (2000). Journal of Geophysical Research: Solid Earth, 105(B5), 11127-11152.
Thivet, S., Carlier, J., Gurioli, L., Di Muro, A., Besson, P., Smietana, M., . . . Nedelec, J.-M. (2022).
Journal of Volcanology and Geothermal Research, 424, 107500.
Thompson, R., & Flower, M. (1971). Earth and Planetary Science Letters, 12(1), 97-107.

6