Modèle d'espace d'état gaussien général en Python

introduction

Le modèle d'espace d'états est une méthode d'analyse de séries chronologiques. Il est hautement interprétable car il peut être modélisé en décomposant les fluctuations au fil du temps en plusieurs composants. Ces dernières années, plusieurs livres ont été publiés qui expliquent le modèle d'espace d'état d'une manière facile à comprendre.

Non expliqué dans cet article

• Qu'est-ce qu'un modèle d'espace d'état et un filtre de Kalman?
• Implémentation R dans KFAS, dlm
• Utilisation de base des modèles de statistiques

Voir les ouvrages de référence et les articles de commentaires ci-dessous.

Que expliquer dans cet article

• Comment implémenter librement un modèle d'espace d'états en utilisant MLEModel of statsmodels
• Application du modèle implémenté à des échantillons de données

Références pour les modèles d'espace d'états

1. [Principes de base de l'analyse des séries temporelles et du modèle d'espace d'état: théorie et mise en œuvre apprises avec R et Stan](https://www.amazon.co.jp/%E6%99%82%E7%B3%BB%E5%] 88% 97% E5% 88% 86% E6% 9E% 90% E3% 81% A8% E7% 8A% B6% E6% 85% 8B% E7% A9% BA% E9% 96% 93% E3% 83% A2% E3% 83% 87% E3% 83% AB% E3% 81% AE% E5% 9F% BA% E7% A4% 8E-R% E3% 81% A8Stan% E3% 81% A7% E5% AD% A6% E3% 81% B6% E7% 90% 86% E8% AB% 96% E3% 81% A8% E5% AE% 9F% E8% A3% 85-% E9% A6% AC% E5% A0% B4 -% E7% 9C% 9F% E5% 93% 89 / dp / 4903814874 / ref = sr_1_1? __ mk_ja_JP =% E3% 82% AB% E3% 82% BF% E3% 82% AB% E3% 83% 8A & mots-clés =% E7% 8A% B6% E6% 85% 8B% E7% A9% BA% E9% 96% 93% E3% 83% A2% E3% 83% 87% E3% 83% AB & qid = 1554681669 & s = passerelle & sr = 8-1) Des trois livres présentés, il contient le moins de formules et est rédigé dans un langage simple. Recommandé pour ceux qui veulent saisir l'image du modèle d'espace d'états pour le moment et l'implémenter dans R. Il est principalement implémenté par KFAS de R, mais il touche également à dlm.

2. [Analyse des séries temporelles à partir du filtre basics-Kalman, MCMC, filtre à particules pratiqué dans R](https://www.amazon.co.jp/%E5%9F%BA%E7%A4%8E%E3%81] % 8B% E3% 82% 89% E3% 82% 8F% E3% 81% 8B% E3% 82% 8B% E6% 99% 82% E7% B3% BB% E5% 88% 97% E5% 88% 86 % E6% 9E% 90-% E2% 80% 95R% E3% 81% A7% E5% AE% 9F% E8% B7% B5% E3% 81% 99% E3% 82% 8B% E3% 82% AB% E3% 83% AB% E3% 83% 9E% E3% 83% B3% E3% 83% 95% E3% 82% A3% E3% 83% AB% E3% 82% BF% E3% 83% BBMCMC% E3% 83% BB% E7% B2% 92% E5% AD% 90% E3% 83% 95% E3% 82% A3% E3% 83% AB% E3% 82% BF% E3% 83% BC-Data-Science- bibliothèque / dp / 4774196460 / ref = pd_sbs_0_1 / 356-7900534-8904633? _encoding = UTF8 & pd_rd_i = 4774196460 & pd_rd_r = 6cba9d2a-5991-11e9-b4dc-811bcb7ea7c6 & pd_rd_w = hpZLR & pd_rd_wg = 8TNR0 & pf_rd_p = ad2ea29d-ea11-483c-9db2-6b5875bb9b73 & pf_rd_r = PQSTK8PGCF3HMAA4XW95 & psc = 1 & refRID = PQSTK8PGCF3HMAA4XW95) Un livre copieux. Il touche également aux filtres Wiener qui ne sont pas mentionnés dans d'autres livres. Vous pouvez également apprendre à implémenter un modèle appliqué avec beaucoup de code. Utilise la bibliothèque dlm de R.

3. [Kalman Filter-Time Series Prediction and State Space Model Using R-](https://www.amazon.co.jp/%E3%82%AB%E3%83%AB%E3%83%9E% E3% 83% B3% E3% 83% 95% E3% 82% A3% E3% 83% AB% E3% 82% BF-% E2% 80% 95R% E3% 82% 92% E4% BD% BF% E3 % 81% A3% E3% 81% 9F% E6% 99% 82% E7% B3% BB% E5% 88% 97% E4% BA% 88% E6% B8% AC% E3% 81% A8% E7% 8A % B6% E6% 85% 8B% E7% A9% BA% E9% 96% 93% E3% 83% A2% E3% 83% 87% E3% 83% AB% E2% 80% 95-% E7% B5% B1% E8% A8% 88% E5% AD% A6One-Point-2 / dp / 4320112539 / ref = pd_sim_0_3 / 356-7900534-8904633? _Encoding = UTF8 & pd_rd_i = 4320112539 & pd_rd_r = 9b8c9003 & pd_rd_r = 9b8c9003 rhkfr & pf_rd_p = b88353e4-7ed3-4da1-bc65-341dfa3a88ce & pf_rd_r = 7Q1J7M8DSYWV6M0ZDF7G & psc = 1 & refRID = 7Q1J7M8DSYWV6M0ZDF7G) Le livre le plus poli avec des explications basées sur des formules mathématiques. Comme son titre l'indique, c'est le plus recommandé pour comprendre le filtre de Kalman. Il y a aussi une explication du filtre à particules. L'exemple de code de KFAS est également inclus.

Tous comprenaient des explications utilisant des formules mathématiques et du code R, ce qui était très utile. D'un autre côté, récemment, une ** implémentation en Python ** peut être requise.

Une bibliothèque appelée ** statsmodels ** est utilisée pour les modèles d'espace d'états en Python, mais il est difficile de comprendre comment l'utiliser car l'explication de la référence officielle est trop petite. Les bases des modèles d'espace d'état utilisant des modèles de statistiques peuvent être trouvées dans cet article.

Cet article décrit une utilisation légèrement plus avancée des modèles de statistiques.

Modèle d'espace d'état gaussien général par modèles de statistiques

Le modèle général d'espace d'état gaussien est exprimé par l'équation suivante. <Mise à jour prévue>

Le modèle de niveau local, qui est le modèle d'espace d'états le plus élémentaire, est un modèle dans lequel toutes les matrices de coefficients sont des matrices unitaires. En spécifiant de manière appropriée la matrice de coefficients, divers modèles tels qu'un modèle de composante de tendance et un modèle de composante périodique peuvent être exprimés. Avec KFAS en R et les modèles de statistiques en Python, vous pouvez créer des modèles en spécifiant directement la matrice de coefficients.

Cependant, il existe un moyen plus simple d'implémenter le modèle de base.

Classe UnobservedComponents

Lors de la mise en œuvre d'un modèle qui incorpore des tendances, des composants périodiques, des composants de régression, etc. dans statsmodels, essentiellement ** Unobserved Components ** Utilisez une classe appelée .structural.UnobservedComponents.html # statsmodels.tsa.statespace.structural.UnobservedComponents).

Cependant, avec UnobservedComponents, le degré de liberté de la modélisation est faible, comme l'incapacité de définir un modèle de composant de régression à l'aide de coefficients variant dans le temps, l'incapacité d'incorporer plusieurs composants périodiques et l'incapacité d'initialiser des états avec des valeurs connues.

Par exemple, l'avantage du filtre de Kalman est qu'il peut faire des ** prédictions séquentielles **, et en pratique il devrait souvent faire des ** prédictions en ligne **. Dans ce cas, le filtre de Kalman est d'abord exécuté sur les données pendant une certaine période de temps, puis les données nouvellement obtenues sont à nouveau filtrées. A ce moment, la valeur initiale de l'état à utiliser pour le filtrage suivant et la dispersion des erreurs de prédiction sont déterminées sur la base de la valeur prédite de l'état au moment final du filtre de Kalman précédent.

Cependant, les composants non observés utilisent par défaut une initialisation distraite approximative (simplement en pensant qu'il n'y a aucun indice sur l'état initial), et ** ne peut pas être initialisé avec des valeurs connues **. Dans ces cas, vous devez utiliser le modèle MLE.

Classe MLEModel

Si vous souhaitez définir votre propre modèle, veuillez consulter ** MLE Model **. https://www.statsmodels.org/stable/generated/statsmodels.tsa.statespace.mlemodel.MLEModel.html?highlight=mlemodels)というクラスを使えば良い。 Il y a Explication sur le site officiel pour savoir comment implémenter le modèle de tendance du système de ligne local. Cependant, ce modèle peut être implémenté avec des composants non observés.

Ici, nous allons introduire l'implémentation d'un modèle qui satisfait les deux points suivants qui ne peuvent pas être réalisés par des composants non observés.

1. Modèle de composante de régression avec coefficient de régression variant dans le temps
2. La méthode d'initialisation peut être sélectionnée

CustomizeMLEmodel.py

import numpy as np
import pandas as pd
import statsmodels.api as sm
import matplotlib.pyplot as plt

"""
Regression Model
"""
class StateSpaceRegression(sm.tsa.statespace.MLEModel):
    def __init__(self, endog, exog, initialization, initial_state=None, initial_state_cov=None):
        # Model order
        # k_state:Dimension d'état
        # k_posdef:Erreur de processus Dimensions de la matrice de dispersion
        k_states = k_posdef = 2

        # Initialize the statespace
        super(StateSpaceRegression, self).__init__(
            endog, exog=exog, k_states=k_states, k_posdef=k_posdef,
            initialization=initialization, initial_state=initial_state, initial_state_cov=initial_state_cov,
            loglikelihood_burn=k_states
        )

        # Initialize the matrices
        #Les variables explicatives sont représentées par une matrice de conception variant dans le temps.
        n_obs = len(exog)
        self.ssm['design'] = np.vstack((np.repeat(1, n_obs), exog))[np.newaxis, :, :]
        
        self.ssm['transition'] = self.ssm['selection'] = np.eye(k_states)

        #auto avec mise à jour.état à ssm_Utilisé pour définir la valeur de cov
        self._state_cov_idx = ('state_cov',) + np.diag_indices(k_posdef)

    @property
    def param_names(self):
        return ['sigma2.measurement', 'sigma2.intercept', 'sigma2.coefficient']

    @property
    #Valeur initiale lors de l'estimation du paramètre le plus probable
    def start_params(self):
        return [np.std(self.endog)]*3

    def transform_params(self, unconstrained):
        return unconstrained**2

    def untransform_params(self, constrained):
        return constrained**0.5

    def update(self, params, *args, **kwargs):
        params = super(StateSpaceRegression, self).update(params, *args, **kwargs)
        
        # Observation covariance
        self.ssm['obs_cov',0,0] = params[0]

        # State covariance
        self.ssm[self._state_cov_idx] = params[1:]

Les points sont les suivants.

1. k_state, k_posdef: Les dimensions du vecteur d'état et de la matrice de dispersion des erreurs de processus. Puisque la dimension de la variable objective (valeur observée) est déduite des données, il n'est pas nécessaire de la donner.
2. mise à jour: utilisé lors de l'estimation du paramètre le plus probable avec la méthode d'ajustement. Prend des paramètres et définit la matrice d'espace d'états appropriée.
3. Matrices d'espaces d'états: Par défaut, toutes les matrices sont des matrices nulles. Notez qu'il ne s'agit pas d'une matrice unitaire. S'il change avec le temps, il doit être défini avec la méthode de mise à jour. Il est stocké dans self.ssm.
4. start params: valeurs initiales pour l'estimation la plus probable des paramètres de distribution
5. initialisation: vous devez définir la moyenne et la variance de la distribution initiale du vecteur d'état. Si la distribution initiale est connue, initialize_known est utilisé, et si elle est inconnue, initialize_approximate_diffuse (initialisation diffuse approximative) est utilisée.
6. transform: Par défaut, l'estimation des paramètres est effectuée sans contraintes, mais en définissant transform_params, elle peut être transformée en paramètres contraints tels que des distributions qui prennent des valeurs positives.
7. param_names: vous pouvez nommer les paramètres que vous souhaitez déduire.

Application aux données réelles

Voici les résultats des expériences utilisant des exemples de données disponibles sur le Site d'assistance dans la troisième référence ci-dessus.

Le premier tracé est la valeur mesurée. $ Y $ (bleu) est prédit par un modèle de composante de régression du coefficient variant dans le temps avec $ x $ (orange) comme variable explicative. Le deuxième graphique représente la composante de niveau et le troisième graphique représente le coefficient de régression. Dans le cas de ces données, on constate que le coefficient de régression tend à augmenter vers la seconde moitié de la période.

Cliquez-ici pour le code.

nintendo = pd.read_csv('data/NINTENDO.csv', index_col=0, parse_dates=[0])
nikkei = pd.read_csv('data/NIKKEI225.csv', index_col=0, parse_dates=[0])

def scaler(s):
    return (s - s.mean()) / s.std()

endog, exog = scaler(nintendo.Close), scaler(nikkei.Close)

#Modèle de coefficient de régression variable dans le temps
mod_reg = StateSpaceRegression(endog, exog, initialization='approximate_diffuse')

res_mod_reg = mod_reg.fit()
res_mod_reg.summary()

#Prévision à un terme de l'état
mu_pred, beta_pred  = res_mod_reg.predicted_state
#État lissé
mu_smoothed, beta_smoothed = res_mod_reg.smoothed_state


fig, axes = plt.subplots(3, 1, figsize=(8, 8), sharex=True)
endog.plot(ax=axes[0], label='y')
exog.plot(ax=axes[0], label='x')
axes[0].set_title('observed variables')
axes[0].legend()

axes[1].plot(endog.index, mu_smoothed)
axes[1].set_title('level component')

axes[2].plot(endog.index, beta_smoothed)
axes[2].set_title('regression coefficient')

Bonus: option de probabilité TensorFlow

Lors de l'utilisation de statsmodels, les paramètres sont estimés ponctuellement par la méthode la plus probable + filtre de Kalman, mais si la modélisation basienne est possible, il est possible de construire et d'estimer un modèle flexible, et la certitude de la prédiction peut être quantifiée naturellement. Il y a. Il semble qu'il soit standard d'utiliser Stan pour implémenter le modèle d'espace d'état bayésien, mais cela nécessite inévitablement un effort de codage. Aussi, bien qu'il y ait du pystan, je pense qu'il est souvent utilisé dans R.

Lors de l'implémentation d'un modèle d'espace d'état bayésien en Python, une nouvelle option est le module sts TensorFlow Probability (TFP). Voir Article officiel.

Cependant, en décembre 2019, la TFP est toujours en cours de développement, elle n'est donc pas encore stable. Attentes pour l'avenir!