示例07 MFDCCA 双序列多分形

"""ex07_mfdcca.py — 第 6 章 MFDCCA 双序列多分形交叉相关

目标：在 DCCA 基础上引入 q 阶，分析两条序列的交叉多分形结构，
得到 h(q)、τ(q)、f(α)。

对应教程：tutorials/06-dcca-mfdcca.md 6.5 节。

运行：
python ex07_mfdcca.py
"""

import numpy as np
import fathon
from fathon import fathonUtils as fu
import matplotlib.pyplot as plt

def main():
np.random.seed(0)
L = 8000
base = np.random.randn(L)
# 两条序列：基序列 + 相关扰动
p1 = fu.toAggregated(base)
p2 = fu.toAggregated(base * 0.7 + np.random.randn(L) * 0.3)

wins = fu.linRangeByStep(10, L // 4, step=5)
q_list = np.arange(-2, 5, 0.5)

mfdcca = fathon.MFDCCA(p1, p2)
n, F = mfdcca.computeFlucVec(wins, q_list, revSeg=True)
h_q, _ = mfdcca.fitFlucVec()
tau = mfdcca.computeMassExponents()
alpha, mfSpect = mfdcca.computeMultifractalSpectrum()

h2 = h_q[np.argmin(np.abs(q_list - 2))]
print(f"MFDCCA h(q=2) = {h2:.4f}")
print(f"Δh = {h_q.max() - h_q.min():.4f}")
print(f"Δα = {alpha.max() - alpha.min():.4f}")

# ── 三联图 ──
fig, axes = plt.subplots(1, 3, figsize=(14, 4.2))

axes[0].plot(q_list, h_q, 'b-o', ms=4)
axes[0].axhline(0.5, color='gray', ls='--', alpha=0.5)
axes[0].set_xlabel('q'); axes[0].set_ylabel('h(q)')
axes[0].set_title('交叉广义 Hurst'); axes[0].grid(alpha=0.3)

axes[1].plot(q_list, tau, 'g-s', ms=4)
axes[1].set_xlabel('q'); axes[1].set_ylabel('τ(q)')
axes[1].set_title('质量指数'); axes[1].grid(alpha=0.3)

axes[2].plot(alpha, mfSpect, 'r-^', ms=4)
axes[2].set_xlabel('α'); axes[2].set_ylabel('f(α)')
axes[2].set_title(f'交叉多分形谱 (Δα={alpha.max()-alpha.min():.3f})')
axes[2].grid(alpha=0.3)

plt.tight_layout()
plt.show()
​

if name == "main":
main()