|
12 | 12 | from dynaris.estimation.diagnostics import acf as compute_acf |
13 | 13 | from dynaris.estimation.diagnostics import standardized_residuals |
14 | 14 | from dynaris.forecast.forecast import ForecastResult, confidence_bands |
15 | | -from dynaris.plotting.style import CMAP, COLORS, create_figure |
| 15 | +from dynaris.plotting.style import COLORS, apply_style, create_figure |
16 | 16 |
|
17 | 17 | # --------------------------------------------------------------------------- |
18 | 18 | # Filtered vs observed |
@@ -359,3 +359,189 @@ def plot_diagnostics( |
359 | 359 | fig.suptitle(title, fontsize=10, fontweight="medium") |
360 | 360 | fig.tight_layout() |
361 | 361 | return fig |
| 362 | + |
| 363 | + |
| 364 | +# --------------------------------------------------------------------------- |
| 365 | +# Panel — combined overview |
| 366 | +# --------------------------------------------------------------------------- |
| 367 | + |
| 368 | + |
| 369 | +def plot_panel( |
| 370 | + filter_result: FilterResult, |
| 371 | + smoother_result: SmootherResult | None, |
| 372 | + forecast_result: ForecastResult | None, |
| 373 | + model: StateSpaceModel, |
| 374 | + component: int = 0, |
| 375 | + level: float = 0.95, |
| 376 | + forecast_levels: tuple[float, ...] = (0.50, 0.80, 0.95), |
| 377 | + n_history: int | None = None, |
| 378 | + n_lags: int = 20, |
| 379 | + title: str = "", |
| 380 | +) -> Figure: |
| 381 | + """Unified panel with filtered, smoothed, forecast, and diagnostics. |
| 382 | +
|
| 383 | + Layout (2 rows x 3 cols): |
| 384 | +
|
| 385 | + - [0, 0] Filtered vs observed |
| 386 | + - [0, 1] Smoothed vs observed |
| 387 | + - [0, 2] Forecast fan chart |
| 388 | + - [1, 0] Standardized residuals |
| 389 | + - [1, 1] Q-Q plot |
| 390 | + - [1, 2] ACF |
| 391 | +
|
| 392 | + Args: |
| 393 | + filter_result: Output of a Kalman filter pass. |
| 394 | + smoother_result: Output of an RTS smoother (optional). |
| 395 | + forecast_result: Output of a forecast (optional). |
| 396 | + model: The state-space model. |
| 397 | + component: Which observation dimension to plot. |
| 398 | + level: Confidence level for filtered/smoothed bands. |
| 399 | + forecast_levels: Confidence levels for forecast fan chart. |
| 400 | + n_history: Historical points to show in forecast panel. |
| 401 | + n_lags: Number of lags for ACF. |
| 402 | + title: Overall figure title. |
| 403 | +
|
| 404 | + Returns: |
| 405 | + The matplotlib Figure. |
| 406 | + """ |
| 407 | + import matplotlib |
| 408 | + import matplotlib.pyplot as plt |
| 409 | + from scipy import stats |
| 410 | + |
| 411 | + fig, axes = plt.subplots(2, 3, figsize=(13, 6.5)) |
| 412 | + fig.set_facecolor("white") |
| 413 | + for ax in axes.flat: |
| 414 | + apply_style(ax) |
| 415 | + |
| 416 | + obs = np.asarray(filter_result.observations[:, component]) |
| 417 | + t = np.arange(len(obs)) |
| 418 | + n_obs = len(obs) |
| 419 | + |
| 420 | + # --- [0, 0] Filtered --- |
| 421 | + ax = axes[0, 0] |
| 422 | + fitted = np.asarray( |
| 423 | + filter_result.filtered_states @ model.F.T |
| 424 | + )[:, component] |
| 425 | + lower_f, upper_f = confidence_bands( |
| 426 | + filter_result.filtered_states @ model.F.T, |
| 427 | + jnp.einsum("ij,tjk,lk->til", model.F, filter_result.filtered_covariances, model.F) |
| 428 | + + model.V[None, :, :], |
| 429 | + level=level, |
| 430 | + ) |
| 431 | + lower_f = np.asarray(lower_f)[:, component] |
| 432 | + upper_f = np.asarray(upper_f)[:, component] |
| 433 | + |
| 434 | + ax.scatter(t, obs, s=4, color=COLORS["observed"], alpha=0.5, zorder=2) |
| 435 | + ax.plot(t, fitted, linewidth=1.0, color=COLORS["secondary"], zorder=3, label="Filtered") |
| 436 | + ax.fill_between(t, lower_f, upper_f, alpha=0.15, color=COLORS["ci_fill"], zorder=1) |
| 437 | + ax.set_title("Filtered", fontsize=8, fontweight="medium") |
| 438 | + ax.legend(fontsize=6, frameon=False) |
| 439 | + |
| 440 | + # --- [0, 1] Smoothed --- |
| 441 | + ax = axes[0, 1] |
| 442 | + if smoother_result is not None: |
| 443 | + smoothed = np.asarray( |
| 444 | + smoother_result.smoothed_states @ model.F.T |
| 445 | + )[:, component] |
| 446 | + lower_s, upper_s = confidence_bands( |
| 447 | + smoother_result.smoothed_states @ model.F.T, |
| 448 | + jnp.einsum( |
| 449 | + "ij,tjk,lk->til", model.F, smoother_result.smoothed_covariances, model.F |
| 450 | + ) |
| 451 | + + model.V[None, :, :], |
| 452 | + level=level, |
| 453 | + ) |
| 454 | + lower_s = np.asarray(lower_s)[:, component] |
| 455 | + upper_s = np.asarray(upper_s)[:, component] |
| 456 | + ax.scatter(t, obs, s=4, color=COLORS["observed"], alpha=0.5, zorder=2) |
| 457 | + ax.plot(t, smoothed, linewidth=1.0, color=COLORS["secondary"], zorder=3, label="Smoothed") |
| 458 | + ax.fill_between(t, lower_s, upper_s, alpha=0.15, color=COLORS["ci_fill_alt"], zorder=1) |
| 459 | + ax.set_title("Smoothed", fontsize=8, fontweight="medium") |
| 460 | + ax.legend(fontsize=6, frameon=False) |
| 461 | + else: |
| 462 | + ax.scatter(t, obs, s=4, color=COLORS["observed"], alpha=0.5, zorder=2) |
| 463 | + ax.plot(t, fitted, linewidth=1.0, color=COLORS["secondary"], zorder=3) |
| 464 | + ax.set_title("Filtered (no smoother)", fontsize=8, fontweight="medium") |
| 465 | + |
| 466 | + # --- [0, 2] Forecast --- |
| 467 | + ax = axes[0, 2] |
| 468 | + if forecast_result is not None: |
| 469 | + fc_mean = np.asarray(forecast_result.mean[:, component]) |
| 470 | + n_fc = len(fc_mean) |
| 471 | + t_fc = np.arange(n_obs, n_obs + n_fc) |
| 472 | + |
| 473 | + hist_obs = obs |
| 474 | + hist_t = t |
| 475 | + hist_fit = fitted |
| 476 | + if n_history is not None: |
| 477 | + hist_obs = obs[-n_history:] |
| 478 | + hist_t = t[-n_history:] |
| 479 | + hist_fit = fitted[-n_history:] |
| 480 | + |
| 481 | + ax.scatter(hist_t, hist_obs, s=4, color=COLORS["observed"], alpha=0.5, zorder=2) |
| 482 | + ax.plot(hist_t, hist_fit, linewidth=0.8, color=COLORS["secondary"], alpha=0.6, zorder=3) |
| 483 | + ax.plot(t_fc, fc_mean, linewidth=1.0, color=COLORS["secondary"], zorder=4, label="Forecast") |
| 484 | + |
| 485 | + cmap = matplotlib.colormaps["Blues"] |
| 486 | + sorted_levels = sorted(forecast_levels, reverse=True) |
| 487 | + n_levels = len(sorted_levels) |
| 488 | + for i, lev in enumerate(sorted_levels): |
| 489 | + lo, hi = confidence_bands( |
| 490 | + forecast_result.mean[:, component], |
| 491 | + forecast_result.covariance[:, component, component], |
| 492 | + level=lev, |
| 493 | + ) |
| 494 | + frac = 0.25 + 0.5 * (i / max(n_levels - 1, 1)) |
| 495 | + ax.fill_between( |
| 496 | + t_fc, np.asarray(lo), np.asarray(hi), |
| 497 | + alpha=0.20, color=cmap(frac), label=f"{int(lev * 100)}%", |
| 498 | + ) |
| 499 | + ax.axvline(n_obs - 0.5, color="#AAAAAA", linewidth=0.4, linestyle="--", zorder=1) |
| 500 | + ax.set_title("Forecast", fontsize=8, fontweight="medium") |
| 501 | + ax.legend(fontsize=6, frameon=False, ncol=2) |
| 502 | + else: |
| 503 | + ax.scatter(t, obs, s=4, color=COLORS["observed"], alpha=0.5) |
| 504 | + ax.set_title("Forecast (none)", fontsize=8, fontweight="medium", color="#AAAAAA") |
| 505 | + |
| 506 | + # --- [1, 0] Standardized residuals --- |
| 507 | + ax = axes[1, 0] |
| 508 | + resids = np.asarray(standardized_residuals(filter_result, model)) |
| 509 | + if resids.ndim > 1: |
| 510 | + resids = resids[:, 0] |
| 511 | + ax.scatter(t, resids, s=3, color=COLORS["secondary"], alpha=0.6) |
| 512 | + ax.axhline(0, color="#888888", linewidth=0.4) |
| 513 | + ax.axhline(2, color=COLORS["ci_fill_alt"], linewidth=0.4, linestyle="--", alpha=0.5) |
| 514 | + ax.axhline(-2, color=COLORS["ci_fill_alt"], linewidth=0.4, linestyle="--", alpha=0.5) |
| 515 | + ax.set_title("Residuals", fontsize=8, fontweight="medium") |
| 516 | + ax.set_ylabel("Std. resid.", fontsize=7) |
| 517 | + |
| 518 | + # --- [1, 1] QQ-plot --- |
| 519 | + ax = axes[1, 1] |
| 520 | + sorted_resids = np.sort(resids) |
| 521 | + n = len(sorted_resids) |
| 522 | + theoretical_q = stats.norm.ppf((np.arange(1, n + 1) - 0.5) / n) |
| 523 | + ax.scatter(theoretical_q, sorted_resids, s=4, color=COLORS["secondary"], alpha=0.6) |
| 524 | + lims = [min(theoretical_q.min(), sorted_resids.min()), |
| 525 | + max(theoretical_q.max(), sorted_resids.max())] |
| 526 | + ax.plot(lims, lims, linewidth=0.7, color=COLORS["ci_fill_alt"], linestyle="--") |
| 527 | + ax.set_xlabel("Theoretical", fontsize=7) |
| 528 | + ax.set_ylabel("Sample", fontsize=7) |
| 529 | + ax.set_title("Q-Q plot", fontsize=8, fontweight="medium") |
| 530 | + |
| 531 | + # --- [1, 2] ACF --- |
| 532 | + ax = axes[1, 2] |
| 533 | + acf_vals = np.asarray(compute_acf(jnp.array(resids), n_lags=n_lags)) |
| 534 | + lags = np.arange(len(acf_vals)) |
| 535 | + ax.bar(lags[1:], acf_vals[1:], width=0.6, color=COLORS["secondary"], alpha=0.6) |
| 536 | + sig = 1.96 / np.sqrt(len(resids)) |
| 537 | + ax.axhline(sig, color=COLORS["ci_fill_alt"], linewidth=0.4, linestyle="--", alpha=0.6) |
| 538 | + ax.axhline(-sig, color=COLORS["ci_fill_alt"], linewidth=0.4, linestyle="--", alpha=0.6) |
| 539 | + ax.axhline(0, color="#888888", linewidth=0.3) |
| 540 | + ax.set_xlabel("Lag", fontsize=7) |
| 541 | + ax.set_ylabel("ACF", fontsize=7) |
| 542 | + ax.set_title("Autocorrelation", fontsize=8, fontweight="medium") |
| 543 | + |
| 544 | + if title: |
| 545 | + fig.suptitle(title, fontsize=11, fontweight="medium", y=1.0) |
| 546 | + fig.tight_layout(h_pad=1.8, w_pad=1.5) |
| 547 | + return fig |
0 commit comments