0. Abstract
- 巨大地震（M>8のもの）は多くの場合、月による潮汐応力による地球の変形が、発生のトリガとなっているらしい（例えば[1]）
- 過去に発生した地震を、統計的に扱った論文もいくつかあった
- しかし、当然といえば当然なのだが、Tidal Phase Angleについて議論していても、Lunar AgeやLunar Phaseを明確に示したものは見つからなかった
- Tidal Phase Angleは、Lunar Ageとほぼ同じような振るまいとなるので、当然といえば当然か
- 一般人が自分でTidal Phase Angleを計算できるとは到底思えない
- しかし、月と太陽の位置なら自分で見ればわかる
- なので、Lunar AgeとEarathquakeの発生状況を可視化した
1. Introduction
このレポートは、過去の地震が統計的にどのような 月の位置とその位相の時に発生したのか を可視化するものである。
月からの潮汐の影響を初めて議論したのは、いまから100年以上も前のSchuster (1987)[1]である。その後、地球潮汐と地震の関係については、特定の震源(群)を対象とした研究としてTanaka (2012)[4]や Ide et al. (2016)[5]といった研究がある。またTsuruoka and Ohtake (1995)[2]は、世界で発生した地震を対象にして地球潮汐との関係を議論したが、非学術領域の一般の人々からみて難解である。そのため、結局のところ月と地震の関係については、世の多くの人から見ると、いまいち要領を得ない。
そこで月や太陽を見れば、その地震が発生する確率を想起できるようにすることが、このレポートの目的である。
2. Data
2.1 Data Source
可視化対象とする地震カタログは、米国地質研究所(USGS: United States Geological Survey)が提供した（と思われる）、Significant Earthquakes, 1965-2016を利用する。このデータは、マグニチュードMが5.5以上のものだけがカタログされている。
USGSだけではなく、日本の気象庁や防災研でも過去の地震カタログを公開している。しかしシンプルな一覧になっておらず、非常に使いにくい。そのため、上述のデータを利用する。
本来、地震はそれぞれ各地点でのプレートテクトニクスによるものなので、地域特性というものが必ず存在する。今回、全地球という広域の地震活動を対象としたため、地域特性を一切考慮していない。そのために、本職の研究者から見たときに何かしら違和感を覚える内容かもしれないことにご留意いただきたい。
2.2 Data Overview
まずはデータを読みこんで、最初の数行を確認する。分析と可視化は、pythonを利用する。
| 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 | import math import datetime import os, sys import numpy as np import pandas as pd DATA_DIR = "/kaggle/input/earthquake-database/" + os.sep # read data file earthquake = pd.read_csv( DATA_DIR+"database.csv", sep=",", parse_dates={'datetime':['Date', 'Time']}, encoding="utf-8", error_bad_lines=False, ) # treating irregular data for idx in [3378,7512,20650]: earthquake.at[idx, "datetime"] = earthquake.at[idx, "datetime"].split(" ")[0] earthquake["datetime"] = pd.to_datetime(earthquake["datetime"], utc=True) earthquake.set_index(["datetime"], inplace=True) earthquake.head() """ Output Latitude Longitude Type Depth Depth Error Depth Seismic Stations Magnitude Magnitude Type Magnitude Error Magnitude Seismic Stations ... Horizontal Error Root Mean Square ID Source Location Source Magnitude Source Status inner moon_dist sun_dist datetime 1965-01-02 13:44:18 19.246 145.616 Earthquake 131.6 NaN NaN 6.0 MW NaN NaN ... NaN NaN ISCGEM860706 ISCGEM ISCGEM ISCGEM Automatic 0.364867 4.128686e+08 1.471089e+11 1965-01-04 11:29:49 1.863 127.352 Earthquake 80.0 NaN NaN 5.8 MW NaN NaN ... NaN NaN ISCGEM860737 ISCGEM ISCGEM ISCGEM Automatic -0.996429 4.065270e+08 1.471063e+11 1965-01-05 18:05:58 -20.579 -173.972 Earthquake 20.0 NaN NaN 6.2 MW NaN NaN ... NaN NaN ISCGEM860762 ISCGEM ISCGEM ISCGEM Automatic 0.947831 4.052391e+08 1.471037e+11 1965-01-08 18:49:43 -59.076 -23.557 Earthquake 15.0 NaN NaN 5.8 MW NaN NaN ... NaN NaN ISCGEM860856 ISCGEM ISCGEM ISCGEM Automatic 0.248578 3.896846e+08 1.471106e+11 1965-01-09 13:32:50 11.938 126.427 Earthquake 15.0 NaN NaN 5.8 MW NaN NaN ... NaN NaN ISCGEM860890 ISCGEM ISCGEM ISCGEM Automatic -0.988605 3.882323e+08 1.471218e+11 """ |
こんな感じのデータとなっていた。
2.3 World Map
地震の震源分布を確認する。
| 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 | import matplotlib.pyplot as plt from mpl_toolkits.basemap import Basemap %matplotlib inline ti = "Map of Earthquake's epicenter duaring 1965-2016" fig = plt.figure(figsize=(18, 18), dpi=96) plt.rcParams["font.size"] = 24 m = Basemap(projection='robin', lat_0=0, lon_0=-170, resolution='c') m.drawcoastlines() m.fillcontinents(color='#606060', zorder = 1) for i in range(5,10,1): #print(i) tmp = earthquake[(earthquake["Magnitude"]>=i)&(earthquake["Magnitude"]<i+1)&(earthquake["Type"]=="Earthquake")] x, y = m(list(tmp.Longitude), list(tmp.Latitude)) points = m.plot(x, y, "o", label=f"Mag.: {i}.x", markersize=0.02*float(i)**3.2, alpha=0.55+0.1*float(i-5)) plt.title(f"{ti}", fontsize=22) plt.legend(bbox_to_anchor=(1.01, 1), loc='upper left', borderaxespad=0, fontsize=18) plt.savefig(OUTPUT_DIR+f'{ti}.png', bbox_inches='tight') |
環太平洋造山帯に、地震が集中していることが明瞭である。っていうか日本ヤバイ。
2.4 Distribution of the Depth
震源の深さ分布を確認する。
| 1 2 3 4 5 6 7 8 9 10 11 12 13 | ti = "Distribution of Earthquake's Depth" plt.figure(figsize=(20, 12), dpi=96) plt.rcParams["font.size"] = 24 for i in range(5,9,1): #print(i) tmp = earthquake[(earthquake["Magnitude"]>=i)&(earthquake["Magnitude"]<i+1)&(earthquake["Type"]=="Earthquake")] plt.hist(tmp["Depth"], bins=50, density=True, histtype='step', linewidth=2.5, label=f"Mag.: {i}.x") plt.legend(bbox_to_anchor=(1.02, 1), loc='upper left', borderaxespad=0, fontsize=18) plt.xlabel("Depth, km") plt.ylabel("Count of Earthquake (Normarized at Total surface=1)") plt.title(f"{ti}") plt.savefig(OUTPUT_DIR+f'{ti}.png', bbox_inches='tight') |
大半の地震は、
- 深さ100㎞以下の表層部分に集中している
- 深さ550-650kmの領域にも、わずかながら発生している
ことがわかる。このグラフだけからはわからないが、深さ550kmを超えるような地震はおそらくプレート境界型ではなく、大陸中央部の大陸地殻内部で発生した地震なのだろう。
潮汐応力や海洋潮汐を起因とする地震の発生震度は、そのメカニズムからしておよそ最大70-80kmまでとなる。以降、可視化の対象とするデータは、80km以浅のものを対象とすることにする。
| 1 2 3 | earthquake = earthquake[earthquake["Depth"]<80] earthquake = earthquake[earthquake["Type"]=="Earthquake"] |
3. Orbit calculation of the Moon and the Sun from the Earth
地球表層部における月からの潮汐応力は、以下の時に大きくなる
- 月―地球間の距離が近い時
- 起潮力は、距離の 3乗に反比例 するため 参考
- 距離の影響が強いため、太陽の影響は月のおよそ45%
- 新月か満月の時
- 太陽と月からの万有引力が同じ方向となり、万有引力のベクトルが揃うため
- 満潮か干潮の前後
- 海面の高さ変化がが、地殻への荷重の変化となるため
そこで、取得した地震リストの各時刻に対して
- 月―地球間の距離
- 太陽―地球間の距離
- 震源地から見た太陽と月の位相角
- 震源地から見た月の方位と高度
を求めることにする。
3.1 Library
天体の位置計算に便利なライブラリであるAstoropyを使用する。
3.2 Eephemeris for Astronomical position calculation
天体の位置計算を行う際の天文暦として、JPL (Jet Propulsion Laboratory:米国ジェット推進研究所)が公開している天文暦のうち、最新のDE432sを使用する。
DE432sは、DE432より短期間（1950-2050）である代わりに、ファイルサイズが小さく（~10MB）、扱いやすい。
3.3 Obtaining of Phase Angle between the Sun and the Moon
月齢を求めたい。
が、月齢そのものが知りたいのではなく、満月と新月のタイミングを知りたい。満月は太陽-地球-月と一直線にならび、新月は太陽-月-地球と一直線に並ぶタイミングとなる。したがって、地球上のある観測点からみた、太陽と月の間の角度をが得られれば、0度の時（太陽と月が同じ方向）は新月、180度の時（太陽と月が反対方向）は満月となることがわかる。
この月と太陽の間の角度を、Phase Angle（位相角）とする。
実際の軌道は3次元的なので、Phase Angleが0度や180度の時に、必ずしも月食や日食になるわけではない。
その計算方法
| 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 | import numpy as np from numpy import linalg as LA import math from astropy.coordinates import EarthLocation, get_body, AltAz, Longitude, Latitude, Angle from astropy.time import Time, TimezoneInfo from astropy import units as u from astropy.coordinates import solar_system_ephemeris solar_system_ephemeris.set('de432s') # 各地震の震央 pos = EarthLocation(Longitude(earthquake["Longitude"], unit="deg"), Latitude(earthquake["Latitude"], unit="deg")) # 位置を求めたい時刻一覧（地震発生時刻のリスト） dts = Time(earthquake.index) # 月の位置（赤道座標から地平座標へ変換） Mpos = get_body("moon", dts).transform_to(AltAz(location=pos)) # 太陽の位置（赤道座標から地平座標へ変換） Spos = get_body("sun", dts).transform_to(AltAz(location=pos)) # phase angle between the Sun and the Moon (rad) SM_angle = Mpos.position_angle(Spos) # 0-360の範囲から0-180の範囲に変換：もともと180-360だった値は180-0へとする earthquake["p_angle"] = [ deg if deg<180 else 360-deg for deg in SM_angle.degree ] # 内積を使用すると、見かけ上、満月や新月の地震が多く見えてしまうので使用しないことにする # Inner Product: converted from Inner Product, 1:New Moon, -1:Full Moon #earthquake["inner"] = [math.cos(rad) for rad in SM_angle.value] earthquake.head() """ outpot Latitude Longitude Type Depth Depth Error Depth Seismic Stations Magnitude Magnitude Type #Magnitude Error Magnitude Seismic Stations ... Root Mean Square ID Source Location Source Magnitude Source Status inner moon_dist sun_dist moon_direction datetime 1965-01-02 13:44:18 19.246 145.616 Earthquake 131.6 NaN NaN 6.0 MW NaN NaN ... NaN ISCGEM860706 ISCGEM ISCGEM ISCGEM Automatic 0.364867 4.128686e+08 1.471089e+11 224.450528 1965-01-04 11:29:49 1.863 127.352 Earthquake 80.0 NaN NaN 5.8 MW NaN NaN ... NaN ISCGEM860737 ISCGEM ISCGEM ISCGEM Automatic -0.996429 4.065270e+08 1.471063e+11 246.918307 1965-01-05 18:05:58 -20.579 -173.972 Earthquake 20.0 NaN NaN 6.2 MW NaN NaN ... NaN ISCGEM860762 ISCGEM ISCGEM ISCGEM Automatic 0.947831 4.052391e+08 1.471037e+11 119.545542 1965-01-08 18:49:43 -59.076 -23.557 Earthquake 15.0 NaN NaN 5.8 MW NaN NaN ... NaN ISCGEM860856 ISCGEM ISCGEM ISCGEM Automatic 0.248578 3.896846e+08 1.471106e+11 346.375786 1965-01-09 13:32:50 11.938 126.427 Earthquake 15.0 NaN NaN 5.8 MW NaN NaN ... NaN ISCGEM860890 ISCGEM ISCGEM ISCGEM Automatic -0.988605 3.882323e+08 1.471218e+11 263.045252 """ |
ここで1つ注意しておきたい。
潮汐応力は月-地球-太陽が一直線上に並ぶ満月や新月のその瞬間が最大となるが、こうした月や太陽の影響が地球表層部に与える影響は、海洋の移動や、剛性体としての地球の挙動が無視できない。
端的に言えば、潮汐応力と潮汐の影響の間にlag(遅れ)が生じる 。このlagは、この満潮や干潮の時刻が場所によって0-8時間ほど遅れる程度のもので、1日以上遅れることはない。
また地球上の月に近い点、そして遠い点の2か所で大潮が発生している。月は27日で地球の周りをまわるが、地球は24時間で自転する。24時間の間に、月の位置はあまり変わらないため、ある地球上の1点において、 満潮と干潮は1日に2回ずつ発生 することも留意。
4. Data Visualization
4.1 Time Series of Earthquakes with its Magnitude
まずは、横軸に時間、縦軸にPhase Angleとした図を示す。
Phase Angleが0度や180度に近い時、すなわち新月や満月の日の周辺で、plotが集中しているように見えなくもない。そこで、Phase Angleのhistogramを描いてみる。
確かに、世界で発生している地震は、0度と180度付近つまり満月と新月の近くで発生するような偏りがありそうだ。ただし、M>=8.0（もしかしたらM>=7.0もか？）のものは、よりはっきりとした傾向があるといえそう。
4.2 Phase Angle and the Distance
4.1の図の横軸を、月―地球間の距離とした図を示す。
距離が0.97と1.05辺りで、plotが集中しているように見える。なので、これをhistogramで確認する。
実際に集中している。要するに、月の公転軌道における近地点と遠地点付近で、集中していることがわかる。
ピークの位置が、潮汐応力が大きくなる近地点付近だけではなく、遠地点でもピークを持つことから、おそらく潮汐応力の変動が大きな時にプレートへの影響が大きくなり、多くの地震を発生させているのだろうと思われる。このことは、Tsuruoka and Ohtake (1995)[2]や、Elizabeth et al. (2014)[3]でも指摘されている。
4.3 Distribution of Moon’s direction
月は約27日で地球を一周するのに対し、地球は24時間で自転する。また、満潮となる位置は月の方向とその真逆の方向の2か所が満月となる。
したがって、もし月からの潮汐応力の影響を強く受けるのであれば、月が南中する時刻あるいは、その12時間前後付近で地震が多く発生することが予想される。
しかし、潮汐応力の変動が大きな時に多くの地震を発生させているのであれば、満月や新月のその瞬間ではなく、数時間の遅れがあるはずである。それはTsuruoka and Ohtake (1995)[2]でも指摘されているように、海洋水の重量がプレートへの荷重として、地震発生トリガとしての影響を及ぼすからである。
場所によるが、地震が多い環太平洋造山帯付近では、大陸が近いこともあり、潮汐応力のphaseより潮位の極大と極小は多くの地点で2-8時間遅れることが観測されている。ただし、太平洋の真ん中にあるハワイでは、ほぼ月が南中するときに満潮を迎える。
地震が発生した時の、月の方位角（北:0deg）の分布を以下に示す。
90度(真東：月の出)と270度(真西：月の入り)付近の2か所にピークがあることがわかる。月が南中した時が、潮汐応力そのものは大きくなるが、その変化が大きくなる位相は90度（6時間）遅れる。さらに、潮汐による海洋荷重の応答も、およそ2-8時間遅れる事実を加味すれば、この図の意味することと矛盾しない。
確かに、地震発生時の月高度の分布をみてみても、0度付近（月の出ないし月の入り）で発生件数が多い。
したがって、地球規模で見たときに、大きな地震は月の出、ないし月の入りの時刻周辺で発生 しているといえる。
ただし、潮汐の遅れは場所ごとに異なるので、本来このような簡便な比較はできないし、Magnitude 8以上の地震では、月高度と地震発生時刻との間に関係性はなさそうだ。
4.4 Seasonal Trend
月別の地震発生件数を確認する。
2011年は、東北地方太平洋沖地震の影響により、3月のMw5.5以上の地震発生件数が突出して多い。そこで、2011年と2011年以外の、月別の1日当たりの地震発生件数を比較した。
| except 2011 | 2011 | |
|---|---|---|
| Jan | 1.145161 | 1.354839 |
| Feb | 1.104218 | 1.322581 |
| Mar | 1.163772 | 7.354839 |
| Apr | 1.183002 | 1.677419 |
| May | 1.176799 | 1.548387 |
| Jun | 1.094293 | 1.387097 |
| Jul | 1.117866 | 2.032258 |
| Aug | 1.209057 | 1.645161 |
| Sep | 1.198511 | 1.290323 |
| Oct | 1.170596 | 1.193548 |
| Nov | 1.197270 | 1.290323 |
| Dec | 1.209057 | 0.903226 |
2011年を除いた年の1日当たりの地震発生回数は、平均1.164、標準偏差0.040となった。全ての月で、2σ の範囲を超えないことから、季節変動はないものと結論付けても良さそうだ。
ただし、この季節変動は非常に強い地域特性があることが報告されている（例えばHeki (2003)[6]）。特定の震源域や断層などに注目すれば、強い季節変動を持つケースもあると思われる。
5. Summary
地震の発生時刻における月と太陽の位置関係を調べた。その結果、地球全体で見たときの大きな地震（Mw>5.5）で以下のようなことがわかった。
- 統計的な観点で、地震の発生トリガとして月（厳密には月軌道の特徴）は地震に強い関連がありそう
- 地球が月と太陽から受ける潮汐応力の変動に対応して、地震の発生件数が増減している
- 年周期（季節変動）は確認できなかった
さらに、以下のような傾向があることも分かった。
- 地震は、月の公転軌道上の近地点ないし遠地点の周辺で発生することが多い
- 地震の発生日は、 満月ないし新月の日の周辺に若干の偏りがあるが、M>=8のような巨大地震は、明らかに新月か満月の周辺が多い
- 地震の発生時刻は、 月の出ないし月の入り周辺の時刻で多い
本来ならば、個々の地震の震源の特徴をもとに、本当に月潮汐に起因した地震かどうかの精査を行う必要がある。例えば、特定の震源を対象とした研究としてTanaka (2012)[4]や Ide et al. (2016)[5]といった研究がある。
が、そこは本職の研究者に託し、ここでは世界の地震に対する統計的な傾向の話にとどめる。
6. Future Prospects
言われるまでもなく、日本は地震大国である。そのため、過去に地震を起因とする津波や火事などを含めて、大きな災害がたびたび発生した。そのため、地震予知に大きな関心を持ち、地震予知連絡会なる組織が政府主導で運営されている。
しかし、地震はその多くがプレート破壊型の地震で、間接的に同じプレートテクトニクスを起因とする断層型地震を含めると、全地震の大半を占める（出典元は失念）。いずれにしても、両者はともに地球プレートの破壊活動で起こるもので、また一般的に破壊活動の将来予測は、そもそも物理学的に不可能である（よほど特殊な条件下なら不可能ではないが、そうした特殊な条件が自然にそろうことはまずない）。
したがって 「地震予知は原理的に不可能である」 が、最近は行政も防災ではなく減災に力点をシフトさせつつある。減災という観点ならば、天気予報のような感覚で発生確率を扱うことができれば、QOLの向上に寄与できる可能性もある（例えば、今週は大地震の確率が高いから、遠出は再来週にしよう、など）。
7. Future Work
ある程度の範囲の特定地域を対象に、同じような可視化をしてみたいと思う。地域によっては、強く月の影響が表れるところもあるかもしれない。
