Note
このページは Jupyter Notebook から生成されました。 ノートブックをダウンロード (.ipynb)
gwexpy TimeSeries 新機能チュートリアル
このノートブックでは、gwpy を拡張した gwexpy の TimeSeries クラスに追加された新しいメソッドとその使い方を解説します。
gwexpy は GWpy との高い互換性を維持しつつ、信号処理、統計解析、他ライブラリとの相互運用性を大幅に強化しています。
目次
環境セットアップ
信号処理と復調 (Hilbert, Phase, Demodulation)
スペクトル解析と相関 (FFT, Transfer Function, xcorr)
ヒルベルト・ファン変換 (HHT)
統計・前処理 (Impute, Standardize, ARIMA, Hurst, Rolling)
リサンプリングと再インデックス (asfreq, resample)
関数によるフィッティング (fit)
相互運用性 (Pandas, Xarray, Torch, and more)
1. 環境セットアップ
まずは必要なライブラリをインポートし、デモ用のサンプルデータを生成します。
[1]:
import matplotlib.pyplot as plt
import numpy as np
from astropy import units as u
from gwexpy.noise.wave import chirp, exponential, gaussian, sine
from gwexpy.plot import Plot
from gwexpy.timeseries import TimeSeries
# SampleData Generate (100Hz done5seconds Data)
fs = 100
duration = 5.0
# Restore 't' for compatibility with downstream cells if they use it
t = np.arange(0, duration, 1 / fs)
# Sensor 1: 10Hz Sine wave + Noise
s1 = sine(duration=duration, sample_rate=fs, frequency=10, amplitude=1.0)
n1 = gaussian(duration=duration, sample_rate=fs, std=0.2)
ts1 = s1 + n1
ts1.name = "Sensor 1"
ts1.override_unit("V")
# Sensor 2: Signal (Frequency time / duration : 5Hz -> 25Hz) + Function Width
s2 = chirp(duration=duration, sample_rate=fs, f0=5, f1=25, t1=duration)
env = exponential(
duration=duration, sample_rate=fs, tau=2.0, decay=False, amplitude=0.2
)
ts2 = s2 * env
ts2.name = "Chirp Signal"
ts2.override_unit("V")
2. 信号処理と復調
gwexpy では、ヒルベルト変換や包絡線、瞬時周波数の計算、さらにロックインアンプのような復調機能が統合されています。
ヒルベルト変換と包絡線
hilbert と envelope を使用します。
[2]:
# minutes Signal (Analytic Signal) Calculate
ts_analytic = ts2.hilbert()
# (Envelope) Calculate
ts_env = ts2.envelope()
plot = Plot(ts2, ts_env, figsize=(10, 4))
ax = plot.gca()
ax.get_lines()[0].set_label("Original (Chirp)")
ax.get_lines()[0].set_alpha(0.5)
ax.get_lines()[1].set_label("Envelope")
ax.get_lines()[1].set_color("red")
ax.get_lines()[1].set_linewidth(2)
ax.legend()
ax.set_title("Hilbert Transform and Envelope")
plt.show()
瞬時位相と瞬時周波数
instantaneous_phase と instantaneous_frequency を使用します。位相のアンラップ (unwrap) や度の単位 (deg) も指定可能です。
[3]:
# Phase (unwrap=True Phase )
phase_rad = ts2.instantaneous_phase(unwrap=True)
phase_deg = ts2.instantaneous_phase(deg=True, unwrap=True)
# Frequency
freq = ts2.instantaneous_frequency()
plot = Plot(phase_rad, freq, separate=True, sharex=True, figsize=(10, 6))
ax = plot.axes
ax[0].set_ylabel("Phase [rad]")
ax[0].set_title("Instantaneous Phase (Unwrapped)")
ax[1].set_ylabel("Frequency [Hz]")
ax[1].set_title("Instantaneous Frequency")
plt.show()
混合と復調 (Mix-down, Baseband, Lock-in)
特定の周波数成分を抽出したり、直流成分へ落とし込むための機能です。
[4]:
# 1. mix_down: Frequency f0 ( Signal )
ts_mixed = ts1.mix_down(f0=10)
# 2. baseband: demodulation + lowpass + resampling BatchRun
ts_base = ts1.baseband(f0=10, lowpass=5, output_rate=20)
# 3. lock_in: (Amplitude Phase )
amp, ph = ts1.lock_in(f0=10, stride=0.1) # 0.1seconds Mean
res_complex = ts1.lock_in(f0=10, stride=0.1, output="complex") # Complex number
print(amp)
plot = Plot(amp, ph, separate=True, sharex=True, figsize=(10, 6))
ax = plot.axes
ax[0].get_lines()[0].set_label("Amplitude")
ax[0].axhline(1.0, color="gray", linestyle="--", label="Theoretical (Sine amp=1)")
ax[0].set_ylabel("Amplitude")
ax[0].legend()
ax[1].get_lines()[0].set_color("orange")
ax[1].get_lines()[0].set_label("Phase [deg]")
ax[1].set_ylabel("Phase [deg]")
plot.figure.suptitle("Lock-in Amplification Result at 10Hz")
plt.show()
TimeSeries([0.97370674, 0.88468624, 0.94466389, 0.91826691,
0.98590416, 0.99147399, 0.81740429, 0.86872436,
0.96546276, 1.00601311, 0.98659591, 1.07083803,
1.11321147, 0.98065663, 0.97148967, 0.97054522,
1.06806132, 0.93631323, 1.09108669, 0.9713755 ,
1.18166385, 0.91119764, 1.09039806, 1.00949188,
1.06136118, 1.22433427, 1.02024191, 1.05812867,
0.96512703, 1.13246248, 0.95906193, 0.95434896,
1.12796849, 0.97984813, 0.97936352, 1.03768126,
1.08577858, 1.09392016, 1.07316951, 1.20109265,
0.94951475, 0.94891906, 0.96796967, 0.93533327,
0.95583234, 0.95080117, 0.91938513, 1.02361956,
0.97574212, 1.03086801],
unit: V,
t0: 0.0 s,
dt: 0.1 s,
name: Sensor 1,
channel: None)
3. スペクトル解析と相関
GWpy の機能を継承しつつ、過渡信号向けの transient モード FFT や、FFT の直接比による伝達関数計算が追加されています。
拡張された FFT
mode="transient" を使用すると、ゼロパディングや高速な長さへの調整、左右個別のパディング指定が可能です。
[5]:
# FFT (GWpy )
fs_gwpy = ts1.fft()
# Transient Mode:
# pad_left/right 。nfft_mode="next_fast_len" Calculate Length 。
fs_trans = ts1.fft(
mode="transient", pad_left=0.5, pad_right=0.5, nfft_mode="next_fast_len"
)
plot = Plot(fs_gwpy.abs(), fs_trans.abs(), yscale="log", xlim=(0, 50), figsize=(10, 4))
ax = plot.gca()
ax.get_lines()[0].set_label("Standard FFT")
ax.get_lines()[1].set_label("Transient FFT (Padded)")
ax.get_lines()[1].set_alpha(0.7)
ax.legend()
ax.set_title("Comparison of FFT modes")
plt.show()
伝達関数 (Transfer Function) と相互相関 (xcorr)
transfer_function では Welch法 (mode="steady") だけでなく、FFT の直接比 (mode="transient") も選択できます。
[6]:
# Transfer Function Calculate (Sensor 1 do Chirp Signal Transfer Function)
tf_welch = ts2.transfer_function(ts1, fftlength=1)
tf_fft = ts2.transfer_function(
ts1, mode="transient"
) # FFTRatio ( Response )
# CrossCorrelation (xcorr)
corr = ts1.xcorr(ts2, maxlag=0.5, normalize="coeff")
plot = Plot(figsize=(10, 6))
ax1 = plot.add_subplot(2, 1, 1)
ax1.semilogy(tf_welch.frequencies, np.abs(tf_welch), label="Welch method")
ax1.semilogy(tf_fft.frequencies, np.abs(tf_fft), label="Direct FFT ratio", alpha=0.5)
ax1.set_title("Transfer Function Magnitude")
ax1.legend()
ax2 = plot.add_subplot(2, 1, 2)
lag = corr.times.value - corr.t0.value
ax2.plot(lag, corr)
ax2.set_title("Cross-Correlation (normalized)")
ax2.set_xlabel("Lag [s]")
plt.show()
STLT (Short-Time Laplace Transform)
STLT (Short-Time Laplace Transform) は、信号の時間変化に伴う局所構造を抽出するための変換です。 gwexpy では STLT の結果は (time × sigma × frequency) の3D変換として TimePlaneTransform / LaplaceGram で表現されます。
以下は、stlt メソッドを使用して STLT を計算し、特定時刻のスライス (Plane2D) を抽出する例です。 この例では frequencies(Hz)を指定し、任意周波数点で STLT を評価します。
[7]:
# Data ( )
import numpy as np
from gwexpy.plot import Plot
t = np.linspace(0, 10, 1000)
data = TimeSeries(np.sin(2 * np.pi * 1 * t), times=t * u.s, unit="V", name="Demo Data")
# STLT Run
# stride: time / duration , window:
freqs = np.array([0.5, 1.0, 1.5]) # Hz
stlt_result = data.stlt(stride="0.5s", window="2s", frequencies=freqs)
print(f"Kind: {stlt_result.kind}")
print(f"Shape: {stlt_result.shape} (Time x Sigma x Frequency)")
print(f"Time Axis: {len(stlt_result.times)} steps")
print(f"Sigma Axis: {len(stlt_result.axis1.index)} bins")
print(f"Frequency Axis: {len(stlt_result.axis2.index)} bins")
# (t=5.0s)
plane_at_5s = stlt_result.at_time(5.0 * u.s)
print(f"Plane at 5.0s shape: {plane_at_5s.shape}")
# Plane2D Verify
print(f"Axis 1: {plane_at_5s.axis1.name}")
print(f"Axis 2: {plane_at_5s.axis2.name}")
Kind: stlt
Shape: (17, 1, 3) (Time x Sigma x Frequency)
Time Axis: 17 steps
Sigma Axis: 1 bins
Frequency Axis: 3 bins
Plane at 5.0s shape: (1, 3)
Axis 1: sigma
Axis 2: frequency
4. Hilbert-Huang Transform (HHT)
非線形・非定常信号解析のためのヒルベルト・ファン変換 (HHT) 機能です。 Empirical Mode Decomposition (EMD) とヒルベルトスペクトル解析を組み合わせます。
[8]:
# HHT (Hilbert-Huang Transform)
# Empirical Mode Decomposition (EMD) Run 、IMFs (Intrinsic Mode Functions)
# : Run PyEMD (EMD-signal) : `pip install EMD-signal`
try:
import os
import warnings
os.environ["TF_CPP_MIN_LOG_LEVEL"] = "2"
warnings.filterwarnings("ignore", category=UserWarning, module="google.protobuf")
warnings.filterwarnings("ignore", message="Protobuf gencode version")
# EMD Run ( Dictionary)
# method="emd" ( EMD) "eemd" (Ensemble EMD)
# ts2 (ChirpSignal )
imfs = ts2.emd(method="emd", max_imf=3)
print(f"Extracted IMFs: {list(imfs.keys())}")
# IMF Plot
sorted_keys = sorted(
[k for k in imfs.keys() if k.startswith("IMF")], key=lambda x: int(x[3:])
)
if "residual" in imfs:
sorted_keys.append("residual")
# gwexpy.plot.Plot BatchPlot
plot_data = [ts2] + [imfs[k] for k in sorted_keys]
plot = Plot(*plot_data, separate=True, sharex=True, figsize=(10, 8))
# Original Set
ax0 = plot.axes[0]
ax0.get_lines()[0].set_label("Original")
ax0.get_lines()[0].set_color("black")
ax0.legend(loc="upper right")
ax0.set_title("Empirical Mode Decomposition")
# IMFs Set
for i, key in enumerate(sorted_keys):
ax = plot.axes[i + 1]
ax.get_lines()[0].set_label(key)
ax.legend(loc="upper right")
plt.show()
except ImportError:
print("EMD-signal not installed. Skipping HHT demo.")
except Exception as e:
print(f"HHT Error: {e}")
Extracted IMFs: ['IMF1', 'IMF2', 'IMF3', 'residual']
5. 統計・前処理
欠損値補完、標準化、ARIMAモデル、ハースト指数、およびローリング統計量が TimeSeries メソッドとして利用可能です。
[9]:
# Data
ts_nan = ts1.copy()
ts_nan.value[100:150] = np.nan
# 1. impute: (interpolation )
ts_imputed = ts_nan.impute(method="interpolate")
# 2. standardize: (z-score, robust )
ts_z = ts1.standardize(method="zscore")
ts_robust = ts1.standardize(method="zscore", robust=True) # Median/IQR
plot = Plot(ts_nan, ts_imputed, figsize=(10, 4))
ax = plot.gca()
ax.get_lines()[0].set_label("with NaNs")
ax.get_lines()[0].set_color("red")
ax.get_lines()[0].set_alpha(0.3)
ax.get_lines()[1].set_label("Imputed")
ax.get_lines()[1].set_linestyle("--")
ax.legend()
ax.set_title("Missing Value Imputation")
plt.show()
ピーク検出 (Find Peaks)
scipy.signal.find_peaks をラップし、時系列データからピークを検出します。
[10]:
# Peak (Find Peaks)
# height, threshold, distance, prominence, width such as
# Unit Threshold
# ts2 (Chirp + Sine) fromPeak
peaks, props = ts2.find_peaks(height=0.0, distance=50)
print(f"Found {len(peaks)} peaks")
if len(peaks) > 0:
print("First 5 peaks:", peaks[:5])
# Plot
plot = ts2.plot(figsize=(12, 4))
ax = plot.gca()
ax.scatter(
peaks.times.value,
peaks.value,
marker="x",
color="red",
s=100,
label="Peaks",
zorder=10,
)
ax.legend(loc="upper right")
ax.set_title("Peak Detection Result")
plt.show()
Found 9 peaks
First 5 peaks: TimeSeries([0.25680508, 0.32974425, 0.44255219, 0.58570943,
0.81881855],
unit: V,
t0: 0.5 s,
dt: 0.5 s,
name: Chirp Signal_peaks,
channel: None)
ARIMAモデルとハースト指数
※ これらの機能には statsmodels, hurst 等のライブラリが必要です。
[11]:
try:
import os
import warnings
os.environ["TF_CPP_MIN_LOG_LEVEL"] = "2"
warnings.filterwarnings("ignore", category=UserWarning, module="google.protobuf")
warnings.filterwarnings("ignore", message="Protobuf gencode version")
# 3. fit_arima: ARIMA(1,0,0)
model = ts1.fit_arima(order=(1, 0, 0))
resid = model.residuals()
forecast, conf = model.forecast(steps=30)
plot = Plot(ts1.tail(100), forecast, figsize=(10, 4))
ax = plot.gca()
ax.get_lines()[0].set_label("Measured")
ax.get_lines()[1].set_label("Forecast")
ax.get_lines()[1].set_color("orange")
# Fill between
ax.fill_between(
conf["lower"].times.value,
conf["lower"].value,
conf["upper"].value,
alpha=0.2,
color="orange",
)
ax.set_title("ARIMA(1,0,0) Fit & Forecast")
ax.legend()
plt.show()
except Exception as e:
print(f"ARIMA skipping: {e}")
try:
import os
import warnings
os.environ["TF_CPP_MIN_LOG_LEVEL"] = "2"
warnings.filterwarnings("ignore", category=UserWarning, module="google.protobuf")
warnings.filterwarnings("ignore", message="Protobuf gencode version")
# 4. hurst / local_hurst: ( Correlation )
h_val = ts1.hurst()
h_detail = ts1.hurst(return_details=True) # Detail
h_local = ts1.local_hurst(window=1.0) # 1seconds
plot = Plot(h_local, figsize=(10, 4))
ax = plot.gca()
ax.get_lines()[0].set_label("Local Hurst")
ax.axhline(h_val, color="red", linestyle="--", label=f"Global H={h_val:.2f}")
ax.set_ylim(0, 0.2)
ax.set_title("Hurst Exponent Analysis")
ax.legend()
plt.show()
except Exception as e:
print(f"Hurst skipping: {e}")
ローリング統計量 (Rolling Statistics)
Pandas と同様の使い勝手で、rolling_mean, std, median, min, max が利用できます。
[12]:
rw = 0.5 * u.s # 0.5seconds
rmean = ts1.rolling_mean(rw)
rstd = ts1.rolling_std(rw)
rmed = ts1.rolling_median(rw)
rmin = ts1.rolling_min(rw)
rmax = ts1.rolling_max(rw)
plot = Plot(ts1, rmean, rmed, figsize=(10, 4))
ax = plot.gca()
ax.get_lines()[0].set_label("Original")
ax.get_lines()[0].set_alpha(0.3)
ax.get_lines()[1].set_label("Rolling Mean")
ax.get_lines()[1].set_color("blue")
ax.get_lines()[2].set_label("Rolling Median")
ax.get_lines()[2].set_color("green")
ax.get_lines()[2].set_linestyle("--")
ax.fill_between(
rmean.times.value,
rmean.value - rstd.value,
rmean.value + rstd.value,
alpha=0.1,
color="blue",
)
ax.legend()
ax.set_title("Rolling Statistics (Window=0.5s)")
plt.show()
6. リサンプリングと再インデックス
asfreq メソッドによる固定グリッドへの割り当てや、resample メソッドでの「時間ビン集計」が可能になりました。
[13]:
# 1. asfreq: Pandas ('50ms' )
ts_reindexed = ts1.asfreq("50ms", method="pad")
Plot(ts1, ts_reindexed)
plt.legend(["Original", ".asfreq('50ms', method='pad')"])
[13]:
<matplotlib.legend.Legend at 0x773893442780>
[14]:
# 2. resample:
# (10Hz) do Signal (GWpy )
ts_sig = ts1.resample(10)
# ('200ms') do 、 time / durationWidth ( )
ts_binned = ts1.resample("200ms") # mean
# Plot doesn't have a direct 'step' method equivalent in args, so we use ax.step or pass to plot
# However, kwarg 'drawstyle'='steps-post' works in plot()? gwpy Plot wraps matplotlib.
# Let's assume we can modify axes after creation or use standard plot for simplicity if appropriate, but step is specific.
# Better strategy: Create Plot instance, then use ax.step
plot = Plot(figsize=(10, 4))
ax = plot.gca()
ax.plot(ts1, alpha=0.6, label="Original (100Hz)")
ax.step(
ts_binned.times,
ts_binned.value,
where="post",
label="Binned Mean (200ms)",
linewidth=2,
)
ax.legend()
ax.set_title("Resampling: Signal Resampling vs Time Binning")
plt.show()
7. 関数によるフィッティング
gwexpy では iminuit をベースとした強力なフィッティング機能を提供しています。GWpy の元のクラスを汚染しないよう、.fit() メソッドはオプトイン方式になっています。
[15]:
from gwexpy.fitting.models import damped_oscillation
# (Function )
# : A=0.5, tau=0.5, f=15, phi=0
result = data.fit(damped_oscillation, A=0.5, tau=0.5, f=15, phi=0)
# Result Display (iminuit )
print(result)
#
# : x_data Datapoints dotime / duration Array
x_data = data.times.value
best_fit = result.model(x_data)
# Plot
plot = data.plot(label="Noisy Signal")
ax = plot.gca()
ax.plot(data.times, best_fit, label="Best Fit", color="red", linestyle="--")
ax.legend()
ax.set_title("Damped Oscillation Fit")
plt.show()
┌─────────────────────────────────────────────────────────────────────────┐
│ Migrad │
├──────────────────────────────────┬──────────────────────────────────────┤
│ FCN = 0.0005703 (χ²/ndof = 0.0) │ Nfcn = 324 │
│ EDM = 0.000176 (Goal: 0.0002) │ │
├──────────────────────────────────┼──────────────────────────────────────┤
│ Valid Minimum │ Below EDM threshold (goal x 10) │
├──────────────────────────────────┼──────────────────────────────────────┤
│ No parameters at limit │ Below call limit │
├──────────────────────────────────┼──────────────────────────────────────┤
│ Hesse ok │ Covariance accurate │
└──────────────────────────────────┴──────────────────────────────────────┘
┌───┬──────┬───────────┬───────────┬────────────┬────────────┬─────────┬─────────┬───────┐
│ │ Name │ Value │ Hesse Err │ Minos Err- │ Minos Err+ │ Limit- │ Limit+ │ Fixed │
├───┼──────┼───────────┼───────────┼────────────┼────────────┼─────────┼─────────┼───────┤
│ 0 │ A │ 1.00 │ 0.06 │ │ │ │ │ │
│ 1 │ tau │ 0 │ 0.06e6 │ │ │ │ │ │
│ 2 │ f │ 1.0000 │ 0.0025 │ │ │ │ │ │
│ 3 │ phi │ -0.00 │ 0.09 │ │ │ │ │ │
└───┴──────┴───────────┴───────────┴────────────┴────────────┴─────────┴─────────┴───────┘
┌─────┬─────────────────────────────────────────────────────┐
│ │ A tau f phi │
├─────┼─────────────────────────────────────────────────────┤
│ A │ 0.00354 -2.4033831e3 3e-6 -0.0001 │
│ tau │ -2.4033831e3 3.78e+09 -8.714e-3 38.287 │
│ f │ 3e-6 -8.714e-3 6.05e-06 -190e-6 │
│ phi │ -0.0001 38.287 -190e-6 0.00796 │
└─────┴─────────────────────────────────────────────────────┘
8. 相互運用性
主要なデータサイエンス・機械学習ライブラリとの相互変換が非常にスムーズです。
[16]:
# Pandas & Xarray
try:
import os
import warnings
os.environ["TF_CPP_MIN_LOG_LEVEL"] = "2"
warnings.filterwarnings("ignore", category=UserWarning, module="google.protobuf")
warnings.filterwarnings("ignore", message="Protobuf gencode version")
df = ts1.to_pandas(index="datetime")
ts_p = TimeSeries.from_pandas(df)
print("Pandas interop OK")
display(df)
except ImportError:
pass
try:
import os
import warnings
os.environ["TF_CPP_MIN_LOG_LEVEL"] = "2"
warnings.filterwarnings("ignore", category=UserWarning, module="google.protobuf")
warnings.filterwarnings("ignore", message="Protobuf gencode version")
xr = ts1.to_xarray()
ts_xr = TimeSeries.from_xarray(xr)
print("Xarray interop OK")
print(xr)
except ImportError:
pass
Pandas interop OK
time_utc
1980-01-06 00:00:19+00:00 0.136702
1980-01-06 00:00:19.010000+00:00 0.916641
1980-01-06 00:00:19.020000+00:00 0.913339
1980-01-06 00:00:19.030000+00:00 0.698927
1980-01-06 00:00:19.040000+00:00 0.307910
...
1980-01-06 00:00:23.950000+00:00 -0.148789
1980-01-06 00:00:23.960000+00:00 -0.602931
1980-01-06 00:00:23.970000+00:00 -1.523759
1980-01-06 00:00:23.980000+00:00 -0.917031
1980-01-06 00:00:23.990000+00:00 -0.661017
Name: Sensor 1, Length: 500, dtype: float64
Xarray interop OK
<xarray.DataArray 'Sensor 1' (time: 500)> Size: 4kB
array([ 0.13670173, 0.9166413 , 0.91333886, ..., -1.5237591 ,
-0.9170307 , -0.66101744], shape=(500,))
Coordinates:
* time (time) datetime64[us] 4kB 1980-01-06T00:00:19 ... 1980-01-06T00:...
Attributes:
unit: V
name: Sensor 1
channel: None
epoch: 0.0
time_coord: datetime
[17]:
# SQLite ( Save)
import sqlite3
with sqlite3.connect(":memory:") as conn:
ts1.to_sqlite(conn, series_id="my_sensor")
ts_sql = TimeSeries.from_sqlite(conn, series_id="my_sensor")
print(f"SQLite interop OK: {ts_sql.name}")
display(conn)
SQLite interop OK: my_sensor
<sqlite3.Connection at 0x773892e114e0>
[18]:
# Deep Learning (Torch)
try:
import os
import warnings
os.environ["TF_CPP_MIN_LOG_LEVEL"] = "2"
warnings.filterwarnings("ignore", category=UserWarning, module="google.protobuf")
warnings.filterwarnings("ignore", message="Protobuf gencode version")
import torch
_ = torch
t_torch = ts1.to_torch()
ts_f_torch = TimeSeries.from_torch(t_torch, t0=ts1.t0, dt=ts1.dt)
print(f"Torch interop OK (Shape: {t_torch.shape})")
display(t_torch)
except ImportError:
pass
Torch interop OK (Shape: torch.Size([500]))
tensor([ 0.1367, 0.9166, 0.9133, 0.6989, 0.3079, -0.0841, -0.4856, -0.9327,
-0.9629, -0.8584, -0.0235, 0.5997, 0.7687, 0.8529, 0.6217, 0.0314,
-0.5548, -0.8344, -1.0038, -0.1451, -0.1633, 0.3889, 1.2103, 0.5661,
0.6103, 0.1592, -0.6597, -0.6484, -1.1797, -0.5220, 0.0160, 0.9038,
0.6236, 0.8553, 0.4024, 0.0122, -0.7108, -0.9686, -0.6863, -0.6871,
0.1749, 0.7761, 0.5115, 1.3691, 0.7149, 0.0465, -0.4601, -1.0190,
-0.8384, -0.3873, -0.3911, 0.5101, 1.0292, 0.9588, 0.4337, 0.2204,
-0.6701, -1.0966, -0.7684, -0.5658, 0.2054, 0.4869, 0.6806, 0.6916,
0.2962, 0.0925, -0.5204, -0.9999, -0.9653, -0.1879, 0.3064, 0.6143,
0.8960, 0.7870, 0.6509, -0.1282, -0.5529, -0.6885, -0.7764, -0.4266,
0.2874, 0.7619, 0.8107, 0.7721, 0.3497, -0.1492, -0.6007, -0.8963,
-1.0040, -0.8043, -0.0340, 0.2960, 1.0429, 1.3018, 0.7348, 0.1753,
-0.4155, -0.9302, -0.8948, -0.3133, 0.1180, 0.4240, 1.3468, 0.9635,
0.3263, -0.1499, -0.6533, -0.6596, -1.0663, -0.4080, 0.1556, 0.5480,
1.1232, 1.0951, 0.3753, -0.1259, -0.5583, -1.1992, -0.9166, -0.5743,
-0.0623, 0.7488, 1.1167, 0.7347, 0.5372, 0.0044, -0.6714, -0.9596,
-1.3882, -0.7172, 0.2482, 0.5632, 0.7518, 1.2144, 0.4542, -0.3655,
-0.7561, -0.7781, -1.0360, -0.3343, 0.1450, 0.7185, 1.0814, 0.7348,
0.1470, 0.2229, -0.5613, -0.8714, -1.2763, -0.3768, 0.2129, 0.7152,
0.8182, 1.1657, 0.4711, -0.2574, -0.3163, -0.9182, -0.9032, -0.5740,
-0.1262, 1.0122, 1.0592, 0.8196, 0.6233, -0.0773, -0.5736, -0.7521,
-0.9501, -1.0797, 0.0291, 0.2220, 0.8995, 0.7233, 0.5450, -0.1013,
-0.5631, -1.0121, -1.2167, -0.4025, -0.2519, 0.7320, 0.9450, 0.7896,
0.6948, 0.0938, -0.6657, -1.1377, -1.0472, -0.8256, -0.0866, 0.5469,
0.8269, 0.8434, 0.4046, 0.1887, -0.3557, -1.0557, -1.1251, -0.6881,
-0.1494, 0.7928, 1.5398, 1.0117, 0.6271, 0.0445, -0.4023, -1.1784,
-0.9176, -0.7092, 0.0366, 0.3882, 1.1468, 0.6869, 0.4952, 0.2551,
-0.8480, -0.7644, -0.7940, -0.5299, -0.0483, 0.7801, 0.8497, 1.0070,
0.7537, -0.0651, -0.7163, -1.1051, -1.0075, -0.6009, -0.1681, 0.6191,
0.7274, 0.9495, 0.6543, -0.1548, -0.6264, -0.9592, -0.9153, -0.9222,
-0.0211, 0.8006, 0.9565, 0.7401, 0.8568, -0.2419, -0.9430, -0.9890,
-1.0288, -0.3337, 0.2888, 0.9933, 1.2823, 1.2736, 0.9331, 0.3349,
-0.5730, -1.2066, -0.7789, -0.5643, 0.0910, 0.5040, 1.0301, 1.0287,
0.3004, -0.0094, -0.8571, -0.9333, -0.9476, -0.6068, -0.2147, 0.5003,
1.0632, 1.2509, 0.6092, 0.0831, -0.7579, -0.5876, -1.0582, -0.6814,
-0.2338, 0.5694, 1.0191, 0.8147, 0.3030, 0.1900, -0.5026, -1.0588,
-1.1938, -0.2229, 0.0154, 0.4889, 0.9706, 1.2214, 0.4746, 0.0732,
-0.7008, -1.1326, -1.0087, -0.9455, 0.2869, 0.4043, 1.0428, 0.9917,
0.0862, -0.1500, -0.4581, -0.8784, -1.1049, -0.6696, 0.0435, 0.6329,
0.5780, 0.7613, 0.4925, -0.2050, -0.7144, -0.9756, -1.1945, -0.5784,
0.0711, 0.5938, 1.1585, 1.1569, 1.0535, 0.0508, -0.3013, -0.7434,
-1.2884, -0.5362, 0.3149, 0.7320, 0.8136, 0.6776, 0.5013, -0.0415,
-0.7761, -0.7534, -1.0752, -0.9316, -0.1174, 0.6804, 1.0546, 0.6063,
0.4915, -0.1157, -0.5058, -1.0389, -1.0827, -0.5210, 0.1408, 0.7555,
1.0042, 0.6564, 0.6033, 0.0163, -0.7819, -1.2529, -0.7168, -0.7656,
-0.1058, 0.5739, 1.5551, 0.7145, 0.4007, 0.0317, -0.6941, -1.1235,
-0.8891, -0.6136, 0.0647, 0.7845, 1.0342, 1.1745, 0.1739, -0.0116,
-0.5755, -0.9884, -1.0122, -0.9535, -0.0319, 0.5532, 0.9905, 1.0591,
0.7648, -0.1131, -0.7542, -1.0242, -0.8643, -0.6847, -0.2936, 0.7235,
1.1228, 1.0940, 0.6417, 0.3456, -0.7260, -1.1559, -1.2386, -0.6276,
-0.1448, 0.6235, 0.5003, 1.0591, 0.4729, -0.2340, -0.8304, -0.8320,
-1.0254, -0.6178, -0.1515, 0.5024, 0.8211, 1.0675, 0.6606, -0.3663,
-0.4069, -0.9435, -0.8441, -0.5528, -0.0453, 0.5804, 1.0866, 0.7835,
0.4905, 0.1617, -0.7779, -0.7664, -0.8323, -0.7719, -0.1674, 0.2714,
0.8693, 0.5408, 0.8183, -0.0442, -0.4756, -0.8131, -1.1349, -0.8260,
-0.2442, 0.6603, 0.9818, 0.7947, 0.6671, -0.2485, -0.2595, -1.0333,
-0.9587, -0.4454, -0.2666, 0.6399, 0.7738, 1.0972, 0.1259, 0.3418,
-0.3092, -1.0487, -0.7799, -0.9347, 0.3438, 0.4086, 0.6605, 1.0398,
0.7587, 0.1601, -0.6990, -0.9220, -0.9434, -0.1776, 0.0176, 0.2776,
0.7415, 1.3080, 0.7845, -0.1660, -0.1796, -1.3537, -0.7655, -0.6591,
0.0401, 0.5844, 0.9035, 1.1136, 0.9327, 0.0809, -0.4139, -0.9517,
-0.5921, -0.5782, 0.0925, 0.7297, 0.8050, 0.7875, 0.1208, -0.1488,
-0.6029, -1.5238, -0.9170, -0.6610], dtype=torch.float64)
[19]:
# Deep Learning (TensorFlow)
try:
import os
import warnings
os.environ.setdefault("TF_CPP_MIN_LOG_LEVEL", "3")
warnings.filterwarnings(
"ignore", category=UserWarning, module=r"google\.protobuf\..*"
)
warnings.filterwarnings(
"ignore", category=UserWarning, message=r"Protobuf gencode version.*"
)
import tensorflow as tf
_ = tf
t_tf = ts1.to_tensorflow()
ts_f_tf = TimeSeries.from_tensorflow(t_tf, t0=ts1.t0, dt=ts1.dt)
print("TensorFlow interop OK")
display(t_tf)
except ImportError:
pass
[20]:
# ObsPy ( )
try:
import os
import warnings
os.environ["TF_CPP_MIN_LOG_LEVEL"] = "2"
warnings.filterwarnings("ignore", category=UserWarning, module="google.protobuf")
warnings.filterwarnings("ignore", message="Protobuf gencode version")
import obspy
_ = obspy
tr = ts1.to_obspy()
ts_f_obspy = TimeSeries.from_obspy(tr)
print(f"ObsPy interop OK: {tr.id}")
display(tr)
except ImportError:
pass
ObsPy interop OK: .Sensor 1..
.Sensor 1.. | 1980-01-06T00:00:00.000000Z - 1980-01-06T00:00:04.990000Z | 100.0 Hz, 500 samples
まとめ
このチュートリアルでは、gwexpy.TimeSeries の主要な拡張機能について学びました:
信号処理: Hilbert変換、復調、拡張FFTモード。
統計: 高度な相関手法 (MIC, Distance Correlation) と ARIMA モデリング。
データクリーニング: 欠損値補完と標準化。
相互運用性: Pandas, Xarray, PyTorch, TensorFlow, ObsPy とのシームレスな変換。
次のステップ
多チャンネルデータ: TimeSeriesMatrix で多数のチャンネルを一括処理する方法を学ぶ。
スペクトル解析: FrequencySeries を詳しく見る。
4Dフィールド: ScalarField で時空間解析を行う。