【数理最適化】結果の可視化と分析

この記事は「【リアルタイム電力予測】需要・価格・電源最適化ダッシュボード構築記」シリーズの十七日目です
Day1はこちら | 全体構成を見る

前回は、数理最適化モデルの全体構造について説明しました。
今回は、その 最適化結果を可視化し、実績と比較しながらモデルの挙動を確認していきます。前回の処理では、以下のような最適化結果ファイル（30分ごとの電源別出力）が出力されています。

コード

可視化コード

import pandas as pd
import numpy as np
import plotly.graph_objects as go
from datetime import datetime
import pytz

JST = pytz.timezone("Asia/Tokyo")

def jst_now_floor_30min() -> datetime:
    """現在時刻をJSTで30分刻みに丸める"""
    now = datetime.now(JST)
    minute = 0 if now.minute < 30 else 30
    return now.replace(minute=minute, second=0, microsecond=0)
now_floor = jst_now_floor_30min()

COLOR_MAP = {
    "lng":           "#F4A6A6",  # パステルレッド（LNG）
    "coal":          "#F7C1C1",  # パステルピンク（石炭）
    "oil":           "#FFD6E8",  # 明るいピンク（石油）
    "th_other":      "#FFEAF2",  # ごく淡いピンク（火力その他）
    "biomass":       "#D9F9B1",  # パステル黄緑（バイオマス）
    "wind_net":      "#B7E4B6",  # パステルグリーン（風力）
    "hydro":         "#B3D9FF",  # パステルブルー（水力）
    "pstorage_gen":  "#CDEEFF",  # 淡い水色（揚水発電）
    "pv_net":        "#FFFACD",  # レモンシフォン（太陽光）
    "battery_dis":   "#FFDAB3",  # パステルオレンジ（蓄電池放電）
    "import":        "#D3D3D3",  # ライトグレー（連系線受電）
    "misc":          "#D8BFD8",  # パステルパープル（その他）
}

def plot_energy_mix(out_df: pd.DataFrame, now_floor):
    df = out_df.copy()
    df["timestamp"] = pd.to_datetime(df["timestamp"])
    df = df.sort_values("timestamp")

    # 数値にしておく（存在する列だけ）
    num_cols = [
        "pv","wind","hydro","coal","oil","lng","th_other", "biomass", "misc"
        "curtail_pv","curtail_wind",
        "pstorage_gen","pstorage_pump",
        "battery_dis","battery_ch",
        "import",
        "reserve","shed","predicted_demand","total_cost"
    ]
    for c in num_cols:
        if c in df:
            df[c] = pd.to_numeric(df[c], errors="coerce").fillna(0.0)

    # --- 純発電（抑制を引いた分） ---
    df["pv_net"]   = np.maximum(df.get("pv", 0)   - df.get("curtail_pv", 0), 0)
    df["wind_net"] = np.maximum(df.get("wind", 0) - df.get("curtail_wind", 0), 0)

    # --- バランス残差チェック ---
    # モデルの需給式：
    #  PV + 風 + 水力 + 火力 + 揚水発電 + 電池放電 + 連系線 + Shed
    #    = 需要 + 揚水揚水 + 電池充電
    gen_total = (
        df.get("pv", 0)
        + df.get("wind", 0)
        + df.get("hydro", 0)
        + df.get("coal", 0)
        + df.get("oil", 0)
        + df.get("lng", 0)
        + df.get("th_other", 0)
        + df.get("biomass", 0)
        + df.get("misc", 0)
        + df.get("pstorage_gen", 0)
        + df.get("battery_dis", 0)
        + df.get("import", 0)
    )
    demand_total = (
        df.get("load_MW", 0)
        + df.get("pstorage_pump", 0)
        + df.get("battery_ch", 0)
    )
    resid = gen_total + df.get("shed", 0) - demand_total
    max_abs_resid = float(np.abs(resid).max())
    if max_abs_resid > 50:  # 50MW超えなら警告
        print(f"[warn] balance residual max |resid| = {max_abs_resid:.1f} MW (確認推奨)")

    fig = go.Figure()

    # --- 積み上げ（供給側：純発電＋揚水発電＋電池放電＋連系線） ---
    stack = "netgen"
    def add_stack(col, name):
        if col in df and df[col].abs().sum() > 0:
            fig.add_trace(go.Scatter(
                x=df["timestamp"], y=df[col],
                name=name, mode="lines",
                stackgroup=stack,
                line=dict(width=0.0),
                fillcolor=COLOR_MAP.get(col, "#cccccc")  # fallback色

            ))

    # 下から積む順序（ベースに近いものを下に）
    add_stack("pstorage_gen",  "揚水(発電)")
    add_stack("hydro",         "水力")
    add_stack("oil",           "石油")
    add_stack("coal",          "石炭")
    add_stack("lng",           "LNG")
    add_stack("pv_net",        "太陽光(純)")
    add_stack("battery_dis",   "蓄電池(放電)")
    add_stack("wind_net",      "風力(純)")
    add_stack("import",        "連系線(流入)")
    add_stack("th_other",      "火力(その他)")
    add_stack("biomass",       "バイオマス")
    add_stack("misc",          "その他")

    # --- 需要ライン（本来の需要） ---
    if "predicted_demand" in df:
        fig.add_trace(go.Scatter(
            x=df["timestamp"], y=df["predicted_demand"], name="需要（予測）",
            mode="lines", line=dict(color="#4FC3F7", width=2.5)
        ))

    # --- 予備力帯（load → load+reserve を塗る） ---
    if "reserve" in df and df["reserve"].abs().sum() > 0:
        # 下側（需要）をダミーで描く
        fig.add_trace(go.Scatter(
            x=df["timestamp"], y=df["predicted_demand"],
            showlegend=False, line=dict(width=0), hoverinfo="skip"
        ))
        # 上側（需要＋予備力）
        fig.add_trace(go.Scatter(
            x=df["timestamp"], y=df["predicted_demand"] + df["reserve"],
            name="予備力帯",
            fill="tonexty", fillcolor="rgba(255,165,0,0.15)",
            line=dict(width=0),
            hovertemplate="%{x|%H:%M}<br>予備力: %{customdata:.0f} MW",
            customdata=df["reserve"]
        ))

    # --- 抑制は別レイヤ（純発電の上に薄く重ねる） ---
    if "curtail_pv" in df and df["curtail_pv"].sum() > 0:
        fig.add_trace(go.Scatter(
            x=df["timestamp"], y=df["pv_net"] + df["curtail_pv"],
            line=dict(width=0), showlegend=False, hoverinfo="skip"
        ))
        fig.add_trace(go.Scatter(
            x=df["timestamp"], y=df["pv_net"],
            name="出力抑制(PV)",
            fill="tonexty", fillcolor="rgba(255,0,0,0.12)",
            line=dict(width=0),
            hovertemplate="%{x|%H:%M}<br>PV抑制: %{customdata:.0f} MW",
            customdata=df["curtail_pv"]
        ))

    if "curtail_wind" in df and df["curtail_wind"].sum() > 0:
        fig.add_trace(go.Scatter(
            x=df["timestamp"], y=df["wind_net"] + df["curtail_wind"],
            line=dict(width=0), showlegend=False, hoverinfo="skip"
        ))
        fig.add_trace(go.Scatter(
            x=df["timestamp"], y=df["wind_net"],
            name="出力抑制(風力)",
            fill="tonexty", fillcolor="rgba(0,0,255,0.12)",
            line=dict(width=0),
            hovertemplate="%{x|%H:%M}<br>風力抑制: %{customdata:.0f} MW",
            customdata=df["curtail_wind"]
        ))

    # --- 供給不足ライン（実際に賄えた需要＝需要 - Shed） ---
    if "shed" in df and df["shed"].max() > 0:
        fig.add_trace(go.Scatter(
            x=df["timestamp"], y=df["predicted_demand"] - df["shed"],
            name="供給不足(実供給)",
            mode="lines",
            line=dict(color="crimson", dash="dot"),
            hovertemplate="%{x|%H:%M}<br>不足: %{customdata:.0f} MW",
            customdata=df["shed"]
        ))

    # --- 総コストは右軸 ---
    if "total_cost" in df:
        fig.add_trace(go.Scatter(
            x=df["timestamp"], y=df["total_cost"],
            name="総コスト",
            mode="lines",
            line=dict(width=2.5, color="#FFCA28", dash="dot"),
            yaxis="y2"
        ))
    
    # 現在時刻の縦線
    fig.add_shape(
        type="line",
        x0=now_floor, x1=now_floor,
        y0=0, y1=1,
        xref="x", yref="paper",
        line=dict(color="red", dash="dot")  # 赤色＋破線
    )

    # 注釈を追加して「Now」と表示
    fig.add_annotation(
        x=now_floor,
        y=1.05,  # yref="paper" の上端に近い位置
        xref="x", yref="paper",
        text="Now",
        showarrow=False,
        font=dict(color="red", size=12),
        align="center"
    )  

    one_hours_ms = 60 * 60 * 1000  # 1時間ごとに縦グリッド
    fig = style_figure(
        fig,
        x_title="時刻",
        y_title="電力 [MW]",
        x_dtick=one_hours_ms,  # グリッドを細かく
        y_dtick=1000,
    )

    # そのあとで、このグラフ固有のレイアウトを上書き
    fig.update_layout(
        # title=title,
        yaxis2=dict(title="コスト", overlaying="y", side="right", showgrid=False),
        hovermode="x unified",
        margin=dict(l=60, r=60, t=60, b=40),
        height=560,
    )

    return fig

スタイル決め

# テーマに合わせた色
PLOT_BG_COLOR = "#0b1633"   # グラフ内の背景（今より少し明るい紺）
PAPER_BG_COLOR = "#020922"  # グラフ外側（アプリ背景と同じでOK）
GRID_COLOR = "#2c3f66"      # グリッド
TEXT_COLOR = "#f5f5f5"      # 文字色

def style_figure(
    fig,
    x_title: str | None = None,
    y_title: str | None = None,
    x_dtick=None,
    y_dtick=None,
):
    """
    Plotly 図に共通のスタイルを当てる。
    - 背景色
    - グリッド色
    - フォントサイズ
    - 軸タイトル
    - 目盛り間隔 (dtick)
    - 凡例
    """
    fig.update_layout(
        paper_bgcolor=PAPER_BG_COLOR,
        plot_bgcolor=PLOT_BG_COLOR,
        font=dict(color=TEXT_COLOR, size=14),  # 全体のデフォルト文字サイズ
        legend=dict(orientation="h", yanchor="bottom", y=1.02, xanchor="left", x=0)
    )

    fig.update_xaxes(
        showgrid=True,
        gridcolor=GRID_COLOR,
        zeroline=False,
        tickfont=dict(size=12),       # x 軸の目盛りの文字
        title_font=dict(size=16),     # x 軸タイトル
    )
    fig.update_yaxes(
        showgrid=True,
        gridcolor=GRID_COLOR,
        zeroline=False,
        tickfont=dict(size=12),       # y 軸の目盛りの文字
        title_font=dict(size=16),     # y 軸タイトル
    )

    if x_title is not None:
        fig.update_xaxes(title_text=x_title)
    if y_title is not None:
        fig.update_yaxes(title_text=y_title)

    # 目盛りの細かさ（任意）
    if x_dtick is not None:
        fig.update_xaxes(dtick=x_dtick)
    if y_dtick is not None:
        fig.update_yaxes(dtick=y_dtick)
    
    return fig

結果

2025年12月14日（天気：雨）

雨の日は太陽光発電がほとんど期待できません。
このような 再エネが弱い日 に、モデルがどのような電源配分を行うかを見ていきます。

数理最適化結果

実績

• 太陽光が弱いため、火力が主に需給を支えている
• 揚水発電が調整弁として動いている
• 全体としては、実績と同程度の供給水準で推移している

考察

LNGが「実績では需要追従、最適化では単調」になっている

実績を見ると、雨の日は LNG 火力が需要の山に合わせて上下しています。
一方、最適化結果では LNG が時間とともになだらかに増加する 形になっています。
この差は、燃料の物理特性というよりも モデルの構造的な仮定 に起因していると考えます。

• モデル側の要因

  • 燃料費が 時間一定の線形コスト

    • ソルバーは「同じ MW をいつ出しても同じコスト」と判断

  • 発電量が 連続変数

    • 起動・停止、最低出力、起動費が存在しない

  • ランプ制約があるため
    → 「できるだけ早めに上げて、そのまま維持」する解が合理的になる

現実では、

• 石炭：止めにくい（最低出力・最低運転時間）
• LNG / 石油：比較的止めやすく、ピーク対応に使われやすい

といった ユニット（台数）単位の制約 が効いていますが、今回のモデルはすべて連続変数のため、その違いが表現できていません。

2025年12月15日（晴れ）

次に、晴れの日の結果を見ていきます。

数理最適化結果

実績

結果

晴れの日にもかかわらず、太陽光発電の山が実績より明らかに低い という問題が見られました(-_-;)

※原因：コードにあらわれていました

当初の設計

1. 過去実績からpv_cap_today（MW）を推定
2. 当日の気象から稼働率cf(t)を計算
3. 発電可能量をpv_avail(t) = pv_cap_today × cf(t)として最適化に渡す

しかし、ここで使っていた 過去の PV 実績自体が、すでに天気の影響を含んでいる という点が問題でした。つまり、「天気を二重に反映した上限」を作ってしまいました。

よって以下のように修正してみました。

1. 過去の PV 実績と当時の全天日射からcapを」逆算（PV実績/cf）
2. 得られた cap の分布から分位点q（例：0.8）を採用
   当日はpv_avail(t) = pv_cap_today × cf_today(t)のみを使う

このように変更した結果、以下のように太陽光発電の山がきちんとできるようになりました。

考察

夕方（PVが落ちる時間）での調整電源

揚水・蓄電池・石油が実績では調整電源としての役割を担います。数理最適化結果を見ると揚水がギザギザな形ではありますが増えています。石油は山なりに戻ってきているので調整電源としての役割を担っています。
蓄電池はそもそも小さいのであまり影響がないようです。

石油が単調になる理由が晴れでも出ているか

雨の日よりかは単調になっていないようです。制約の範囲内で柔軟に動いてくれています。

誤差の分析

晴れの日の誤差を分析します。

時系列比較（最適化：実線 / 実績：点線）

いくつかの電源をpick upして描画します。

太陽光発電

→ 全天日射の予報に沿った形状を概ね再現できています

風力発電

→ こちらは短周期の変動は捉えきれていません

水力発電

→ モデルでは昼に 0 近くまで落ちるが、実績では緩やかに推移しています

LNG火力

→ 帯グラフで可視化した際はそこまで気になりませんでしたが、朝夕に山がある実績に対し、最適化は単調増加になっています

誤差ヒートマップ（opt − actual）

• LNG：昼夜を通して正の偏り
• Import：時間帯で符号が変化

電源別誤差指標

以下3つを可視化しています。

• MAE：平均的なズレの大きさ
• RMSE：大外しを強調した誤差
• Bias：出しすぎ／出さなすぎの傾向

• LNG / Import の Bias が大きい
  → 調整をこれらに担わせる設計になっている
• PV は Bias が小さい
  → 山の形状は再現できている

需給バランスの比較

青色が数理最適化で表した供給曲線、赤点線が予測した電力需要、緑破線が実際の電力需要を表しています。

昼間に供給が需要（予測）を上回る場面があり、これは PV 上限制約・揚水・予備力の組み合わせの影響と考えられます。

明日

今日で数理最適化の章は終わりです。将来の気象や需要の予測を入れた数理最適化はなかなか難しかったです。なんとか実績値と比べても遜色ない程度にまで仕上げられてよかったです。

次回からは、この結果を 毎日自動更新するダッシュボード を構築するためのデータ取得・パイプライン設計について書いていきます！