3. ünite / 8

Kalman süzgecine giriş: modele mi ölçüme mi güvenmeli

Tahmin ve güncelleme adımları, Q ve R'nin anlamı, kazancın oturması ve tek boyutta Kalman ile üstel ortalamanın akrabalığı.

Bu ünitenin sonunda

Tek boyutlu Kalman süzgecinin tahmin ve güncelleme adımlarını yazmak
Q ve R'nin süzgecin davranışını nasıl değiştirdiğini açıklamak
Durağan kazancı hesaplamak ve aynı kazançlı üstel ortalamayla karşılaştırmak
Kalman'ın ne zaman EMA'dan fazlasını verdiğini söylemek

Ders notları

Bir sensörü okurken iki bilgi kaynağın vardır: model ("sıcaklık on dakikada pek değişmez") ve ölçüm ("DHT22 şu an 24,3 diyor"). İkisi de hatalıdır. Kalman süzgeci her adımda ikisini, belirsizlikleriyle orantılı olarak karıştırır.

Model ve iki gürültü

x_k = x_k−1 + w, Var(w) = Q z_k = x_k + v, Var(v) = R

Q (süreç gürültüsü): gerçek değer bir adımda ne kadar değişebilir? Sera yavaş ısınır, Q küçüktür.
R (ölçüm gürültüsü): sensöre ne kadar güvenebilirsin? DHT22 için yaklaşık 0,25² = 0,0625 °C².

İki adım

Tahmin: x aynı kalır, belirsizlik büyür: P ← P + Q.
Güncelleme: kazanç K = P / (P + R). Tahmin yeniliğin K katı kadar ölçüme kayar: x ← x + K(z − x). Belirsizlik küçülür: P ← (1 − K)P.

K 0'a yakınsa süzgeç modele güvenir (yumuşak ama geç tepki), 1'e yakınsa ölçüme güvenir (hızlı ama gürültülü). Ölçüm yoksa (NaN) yalnız tahmin adımı yapılır: P büyür, sonraki ölçüme daha çok güvenilir. Bu, eksik veriyi doldurmanın ilkeli bir yoludur.

Kazanç oturur

Q ve R sabitse P birkaç adımda sabit bir değere, K de durağan kazanca oturur:

P_∞ = (Q + √(Q² + 4QR)) / 2 K_∞ = P_∞ / (P_∞ + R)

Q = 0,02 ve R = 0,0625 için K_∞ ≈ 0,43'tür. Oturduktan sonra Kalman güncellemesi x += K(z − x) tam olarak α = K olan bir üstel ortalamadır. Tek boyutta, sabit gürültüde Kalman'ın sana kattığı şey α'yı tahminle değil, fiziksel olarak anlamlı iki sayıdan (sensör belgesindeki gürültü ve sürecin hızı) türetmektir. Sera verisinde bu süzgeç gerçek sıcaklığa karşı yaklaşık 0,31 °C RMSE verir.

Kalman ne zaman fazlasını verir?

Değişken R: sensör bazen daha gürültülüyse (ör. GPS uydu sayısı düştüğünde) K her adımda kendini ayarlar.
Eksik ölçüm: P büyüyerek belirsizliği taşır.
Çok boyut ve füzyon: konum ve hız birlikte, ya da barometre ile ivmeölçer birlikte. Drone dersindeki füzyon bunun genelleşmesidir.
Yenilik izleme: z − x yeniliği beklenen dağılımın dışına çıkarsa model bozulmuştur. 5. günkü ısıtıcı arızasının ilk üç saatinde ortalama yenilik +0,68 °C'dir, bir gün önce aynı saatlerde −0,17 °C: bu, anomali uyarısının tohumudur.

Simülasyon

1B Kalman süzgeci: Q ve R ile oyna

Gürültülü bir ölçüm, bir model ve iki sayı: süreç gürültüsü Q (değer ne kadar hızlı değişebilir) ve ölçüm gürültüsü R (sensöre ne kadar güvenilir). Kazanç K'nın nasıl oturduğunu ve aynı kazançlı üstel ortalamayla farkı gör.

K∞ = 0,426 (aynı değerli EMA α'sı) Gerçeğe göre RMSE: ham 0,490 · Kalman 0,339 · EMA(α = K∞) 0,344 °C

Süreç gürültüsü Q2.0e-2Ölçüm gürültüsü R6.3e-2Sensör gürültüsü ×2,0

Sensör σ ≈ 0,50 °C → uygun R = σ² ≈ 0,250

Kod

Tek boyutlu Kalman süzgeci (C++ yapı)

// Model: x_k = x_(k-1) + w, Var(w) = q;  ölçüm z_k = x_k + v, Var(v) = r
struct Kalman1D {
  float q, r;          // süreç ve ölçüm gürültüsü (°C²)
  float x = NAN, p = 1.0f, k = 0.0f;

  float update(float z) {
    if (isnan(x)) { x = z; return x; }   // ilk ölçümle başla
    p += q;                              // tahmin: belirsizlik büyür
    if (isnan(z)) { k = 0; return x; }   // ölçüm yoksa yalnız tahmin
    k = p / (p + r);                     // kazanç: modele mi ölçüme mi güven
    x += k * (z - x);                    // güncelle: yeniliğin K katı
    p *= (1.0f - k);
    return x;
  }
};

// DHT22 için: r = 0,25² (belgedeki gürültü), q = 0,02 (sera yavaş ısınır)
Kalman1D kf{0.02f, 0.0625f};

Python: durağan kazanç ve EMA eşdeğerliği

import math

def steady_gain(q, r):
    p = (q + math.sqrt(q * q + 4 * q * r)) / 2   # tahmin öncesi P sonsuzda
    return p / (p + r)

for q in (0.002, 0.02, 0.2):
    K = steady_gain(q, 0.0625)
    print("q=%.3f  K=%.3f  ~ EMA alpha" % (q, K))  # q büyüdükçe K büyür: daha az yumuşatma

Kodu tarayıcıda çalıştırPython · NumPy · matplotlib · Pyodide

Python ve NumPy bu sayfada, tarayıcında çalışır (Pyodide). Kodu değiştirip yeniden çalıştırabilirsin; dersin veri setleri çalışma klasöründe hazır. Kod matplotlib kullanırsa grafikler çıktının altında görünür, yazdığı dosyalar (ör. egitim_gecmisi.csv) indirilebilir.

kalman.pypython

import numpy as np
D = np.genfromtxt('sera-sensorleri.csv', delimiter=',', names=True)  # 7 gün, 10 dakikada bir
t = D["sicaklik"].copy()
t[(t > 60) | (t < -20)] = np.nan              # imkânsız değerleri at
def kalman(z, q, r):
    x, p = z[~np.isnan(z)][0], 1.0
    est, gain = [], []
    for v in z:
        p += q                                # tahmin
        k = 0.0
        if not np.isnan(v):
            k = p / (p + r); x += k * (v - x); p *= 1 - k  # güncelle
        est.append(x); gain.append(k)
    return np.array(est), np.array(gain)
r = 0.25 ** 2                                 # DHT22 gürültüsü ~0,25 °C
for q in (0.002, 0.02, 0.2):
    est, K = kalman(t, q, r)
    Ks = (q + np.sqrt(q * q + 4 * q * r)) / 2
    Ks = Ks / (Ks + r)                        # durağan kazanç formülü
    print("q=%.3f  K(ilk 5)=%s  K(son)=%.3f  formul=%.3f" % (q, np.round(K[:5], 3), K[K > 0][-1], Ks))  # NaN ölçümde K = 0
est, K = kalman(t, 0.02, r)
a = K[K > 0][-1]; y = np.empty(len(t)); y[0] = est[0]
for k in range(1, len(t)):                    # aynı kazançlı EMA
    y[k] = y[k - 1] if np.isnan(t[k]) else y[k - 1] + a * (t[k] - y[k - 1])
print("Kalman ile EMA arasindaki en buyuk fark (ilk 20 sonrasi): %.4f C" % np.max(np.abs(est - y)[20:]))
inn = t - np.roll(est, 1)                     # yenilik: ölçüm eksi önceki tahmin
gun = (D["dakika"] // 1440).astype(int); h = D["saat"]
ramp = (D["ariza"] == 1) & (h < 16)           # 5. gün 13:00-16:00, ısıtıcı takılı
normal = (gun == 3) & (h >= 13) & (h < 16)    # 4. gün aynı saatler
print("ortalama yenilik  ariza: %+.2f C   normal gun: %+.2f C" % (np.nanmean(inn[ramp]), np.nanmean(inn[normal])))  # model "sabit" der, ısıtıcı ısıtır

Kendi CSV’ni yükle

Çıktı

Çıktı burada görünecek.

Kendin dene

1B Kalman simülasyonunda R'yi 0,0625'te tut, Q'yu 0,002'den 0,2'ye çıkar: tahmin eğrisi ve kazanç K nasıl değişiyor? Sonra Q'yu sabit tutup R'yi 10 kat büyüt ve aynı etkiyi ters yönden elde ettiğini göster. Kazancın oturduğu değeri not et ve aynı α ile çizilen üstel ortalamanın Kalman eğrisinin üstüne bindiğini doğrula. Python kodunda K(son) ile formül değerini karşılaştır ve arızanın ilk saatlerindeki ortalama yeniliği normal günle karşılaştır.

Sık yapılan hatalar

Q ve R'yi gerçek birimlerle düşünmemek: R sensörün varyansıdır (°C²), standart sapması değil.
Kalman'ı aykırı değer temizliğinin yerine koymak: 85 °C'lik sıçrama, K kadar ağırlıkla tahmini bozar; önce aralık ya da Hampel denetimi.
Başlangıç P'sini çok küçük seçmek: süzgeç ilk ölçümlere güvenmez ve yavaş yakınsar.
Model yanlışken Q'yu küçültmeye devam etmek: süzgeç gerçek değişimi gürültü sanar ve geride kalır.

Yazılım mühendisliği bağlantısı

Kalman süzgeci küçük bir durum makinesidir: durum (x, P), girdi (z) ve iki saf geçiş (tahmin, güncelleme).

Parametreleri (Q, R, P₀) koda gömme; yapılandırma olarak ver ve neden o değeri seçtiğini yorumla.
Özellik testi: sabit girdi verildiğinde tahmin o değere yakınsar; R → ∞ iken K → 0; Q → ∞ iken K → 1.
Durağan kazanç formülünü uzun bir döngünün son K değeriyle karşılaştıran bir test, iki kodun da doğru olduğunu kanıtlar.

Kendini sına

R (ölçüm gürültüsü) büyütülürse kazanç K ne olur?
Ölçüm NaN geldiğinde 1B Kalman ne yapar?
Q ve R sabitken, kazanç oturduktan sonra 1B Kalman neye eşdeğerdir?
Isıtıcı arızasının ilk saatlerinde Kalman yeniliği (z − x) nasıl davranır?