Estimation des paramètres avec le filtre de Kalman

Dans les prévisions météorologiques, l'estimation des paramètres est effectuée à l'aide du filtre de Kalman pour convertir la sortie du modèle en variables de prévision (par exemple, la température maximale). La méthode est décrite dans l'Agence météorologique (1997).

Estimation des paramètres

Quand il n'y a pas de dépendance au temps

Ici, estimons $ a $ et $ b $ dans l'équation linéaire $ y = ax + b $. Utilisez la formule $ y_i = a x_i + b + v_i $ pour générer des données simulées et les examiner.

Soit $ a = 2 $, $ b = 5 $, $ x $ un nombre aléatoire uniforme de -5 à 5, et $ v_i $ donne un nombre aléatoire gaussien avec une moyenne de 0 et un écart type de 2.

KLM.py



def update(P, C, R, x_hat, obs, I):
  """
  P:Matrice de covariance d'erreur
  C:Matrice des numéros du système d'observation
  R:Matrice de dispersion du bruit observé
  """
  #Gain de Kalman
  G = P * C / (C.T * P * C + R)
  x_hat = x_hat + G * (obs - C.T * x_hat)
  P = (I -  G * C.T) * P
  return x_hat, P

a = 2.
b = 5.

x = np.random.uniform(-5, 5, 365)
v = np.random.normal(0, 2, 365)
y = []
for i in range(365):
    y.append(a * x[i] + b + v[i])
    
A = np.mat([1])
P = np.mat([[1, 0], [0, 1]])
R = np.mat([1])
I = np.identity(2)
x_hat = np.mat([[0], [0]])

for i in range(365):
  C = np.mat([[x[i]], [1]])
  obs = np.mat([y[i]])
  x_hat, P = update(P, C, R, x_hat, obs, I)

Pour illustrer les estimés $ a $ et $ b $, KLM.png On peut voir que $ a $ et $ b $ sont estimés sur une courte période de temps.

S'il y a une dépendance au temps

Lorsque le modèle change en cours de route ou lorsque les paramètres varient selon les saisons. Cette fois, lorsque 2 est ajouté à chaque coefficient le 180e jour.

KLM2.py


import numpy as np

def update(P, C, R, x_hat, obs, I):
  """
  P:Matrice de covariance d'erreur
  C:Matrice des numéros du système d'observation
  R:Matrice de dispersion du bruit observé
  """
  #Gain de Kalman
  G = P * C / (C.T * P * C + R)
  x_hat = x_hat + G * (obs - C.T * x_hat)
  P = (I -  G * C.T) * P
  return x_hat, P

a = 2.
b = 5.

x = np.random.uniform(-5, 5, 365)
v = np.random.normal(0, 2, 365)
y = []
for i in range(365):
    if i < 180:
        y.append(a * x[i] + b + v[i])
    else:
        y.append((a + 2) * x[i] + (b + 2) + v[i])

A = np.mat([1])
P = np.mat([[1, 0], [0, 1]])
R = np.mat([1])
S = 0.1
I = np.identity(2)
x_hat = np.mat([[0], [0]])

for i in range(365):
  C = np.mat([[x[i]], [1]])
  obs = np.mat([y[i]])
  P = P + S ** 2
  x_hat, P = update(P, C, R, x_hat, obs, I)

KLM2.png

Il a suivi le changement avec brio.

Les références

Agence météorologique (1997) Dernier guide des prévisions numériques Chapitre 5 Système d'application 1 Filtre de Kalman pp.66 --78

Recommended Posts

Estimation des paramètres avec le filtre de Kalman
Détection de point de changement avec filtre de Kalman
Filtre de Kalman non linéaire (1)
Réglage des paramètres avec luigi (2)
Réglage des paramètres avec luigi
Créer un filtre avec scipy
Utiliser DATE_FORMAT avec le filtre SQLAlchemy
[Python] Filtrer les feuilles de calcul avec gspread
Filtre de Kalman que vous pouvez comprendre
Implémentation de l'estimation des paramètres HMM en python
Créer un filtre avec un modèle django
Estimation des paramètres moteurs à partir des points correspondants