diff --git a/doc/source/conf.py b/doc/source/conf.py index e37f9990..81df660c 100644 --- a/doc/source/conf.py +++ b/doc/source/conf.py @@ -18,7 +18,8 @@ # import os import sys -sys.path.insert(0, os.path.abspath('../..')) + +sys.path.insert(0, os.path.abspath("../..")) # -- General configuration ------------------------------------------------ @@ -29,54 +30,54 @@ html_theme = "sphinx_rtd_theme" html_theme_path = [sphinx_rtd_theme.get_html_theme_path()] -#collapse_navigation = True +# collapse_navigation = True # Add any Sphinx extension module names here, as strings. They can be # extensions coming with Sphinx (named 'sphinx.ext.*') or your custom # ones. extensions = [ - 'sphinx.ext.autodoc', - 'sphinx.ext.mathjax', - 'sphinx.ext.githubpages', - 'sphinx.ext.napoleon', + "sphinx.ext.autodoc", + "sphinx.ext.mathjax", + "sphinx.ext.githubpages", + "sphinx.ext.napoleon", ] # Add any paths that contain templates here, relative to this directory. -templates_path = ['_templates'] +templates_path = ["_templates"] # The suffix(es) of source filenames. # You can specify multiple suffix as a list of string: # # source_suffix = ['.rst', '.md'] -source_suffix = '.rst' +source_suffix = ".rst" # The encoding of source files. # # source_encoding = 'utf-8-sig' # The master toctree document. -master_doc = 'index' +master_doc = "index" # General information about the project. -project = u'PyDLM' -copyright = u'2024, Xiangyu Wang' -author = u'Xiangyu Wang' +project = "PyDLM" +copyright = "2024, Xiangyu Wang" +author = "Xiangyu Wang" # The version info for the project you're documenting, acts as replacement for # |version| and |release|, also used in various other places throughout the # built documents. # # The short X.Y version. -version = u'0.1.1.13' +version = "0.1.1.13" # The full version, including alpha/beta/rc tags. -release = u'0.1.1.13' +release = "0.1.1.13" # The language for content autogenerated by Sphinx. Refer to documentation # for a list of supported languages. # # This is also used if you do content translation via gettext catalogs. # Usually you set "language" from the command line for these cases. -language = 'en' +language = "en" # There are two options for replacing |today|: either, you set today to some # non-false value, then it is used: @@ -112,7 +113,7 @@ # show_authors = False # The name of the Pygments (syntax highlighting) style to use. -pygments_style = 'sphinx' +pygments_style = "sphinx" # A list of ignored prefixes for module index sorting. # modindex_common_prefix = [] @@ -163,7 +164,7 @@ # Add any paths that contain custom static files (such as style sheets) here, # relative to this directory. They are copied after the builtin static files, # so a file named "default.css" will overwrite the builtin "default.css". -html_static_path = ['_static'] +html_static_path = ["_static"] # Add any extra paths that contain custom files (such as robots.txt or # .htaccess) here, relative to this directory. These files are copied @@ -243,34 +244,30 @@ # html_search_scorer = 'scorer.js' # Output file base name for HTML help builder. -htmlhelp_basename = 'PyDLMdoc' +htmlhelp_basename = "PyDLMdoc" # -- Options for LaTeX output --------------------------------------------- latex_elements = { - # The paper size ('letterpaper' or 'a4paper'). - # - # 'papersize': 'letterpaper', - - # The font size ('10pt', '11pt' or '12pt'). - # - # 'pointsize': '10pt', - - # Additional stuff for the LaTeX preamble. - # - # 'preamble': '', - - # Latex figure (float) alignment - # - # 'figure_align': 'htbp', + # The paper size ('letterpaper' or 'a4paper'). + # + # 'papersize': 'letterpaper', + # The font size ('10pt', '11pt' or '12pt'). + # + # 'pointsize': '10pt', + # Additional stuff for the LaTeX preamble. + # + # 'preamble': '', + # Latex figure (float) alignment + # + # 'figure_align': 'htbp', } # Grouping the document tree into LaTeX files. List of tuples # (source start file, target name, title, # author, documentclass [howto, manual, or own class]). latex_documents = [ - (master_doc, 'PyDLM.tex', u'PyDLM Documentation', - u'Xiangyu Wang', 'manual'), + (master_doc, "PyDLM.tex", "PyDLM Documentation", "Xiangyu Wang", "manual"), ] # The name of an image file (relative to this directory) to place at the top of @@ -310,10 +307,7 @@ # One entry per manual page. List of tuples # (source start file, name, description, authors, manual section). -man_pages = [ - (master_doc, 'pydlm', u'PyDLM Documentation', - [author], 1) -] +man_pages = [(master_doc, "pydlm", "PyDLM Documentation", [author], 1)] # If true, show URL addresses after external links. # @@ -326,9 +320,15 @@ # (source start file, target name, title, author, # dir menu entry, description, category) texinfo_documents = [ - (master_doc, 'PyDLM', u'PyDLM Documentation', - author, 'PyDLM', 'One line description of project.', - 'Miscellaneous'), + ( + master_doc, + "PyDLM", + "PyDLM Documentation", + author, + "PyDLM", + "One line description of project.", + "Miscellaneous", + ), ] # Documents to append as an appendix to all manuals. diff --git a/examples/unemployment_insurance/main.py b/examples/unemployment_insurance/main.py index 9b1b4328..da249217 100644 --- a/examples/unemployment_insurance/main.py +++ b/examples/unemployment_insurance/main.py @@ -1,4 +1,4 @@ -#================================================================== +# ================================================================== # # This example use the data from Google data science blog post # http://www.unofficialgoogledatascience.com/2017/07/fitting-bayesian-structural-time-series.html @@ -6,57 +6,60 @@ # Make sure your working directly is pointing to the folder where # `main.py` locates and `data.csv` is in the same folder. # -#================================================================== +# ================================================================== # Read data file import os + this_dir = os.getcwd() DATA_PATH = os.path.join(this_dir, "data.csv") -data_file = open(DATA_PATH, 'r') +data_file = open(DATA_PATH, "r") -variables = data_file.readline().strip().split(',') +variables = data_file.readline().strip().split(",") data_map = {} for var in variables: data_map[var] = [] for line in data_file: - for i, data_piece in enumerate(line.strip().split(',')): + for i, data_piece in enumerate(line.strip().split(",")): data_map[variables[i]].append(float(data_piece)) # Extract and store the data. time_series = data_map[variables[0]] -features = [[data_map[variables[j]][i] for j in range(1, len(variables)) ] - for i in range(len(time_series))] +features = [ + [data_map[variables[j]][i] for j in range(1, len(variables))] + for i in range(len(time_series)) +] # Plot the raw data import matplotlib.pyplot as plt import pydlm.plot.dlmPlot as dlmPlot -dlmPlot.plotData(range(len(time_series)), - time_series, - showDataPoint=False, - label='raw_data') -plt.legend(loc='best', shadow=True) + +dlmPlot.plotData( + range(len(time_series)), time_series, showDataPoint=False, label="raw_data" +) +plt.legend(loc="best", shadow=True) plt.show() # Build a simple model from pydlm import dlm, trend, seasonality # A linear trend -linear_trend = trend(degree=1, discount=0.95, name='linear_trend', w=10) +linear_trend = trend(degree=1, discount=0.95, name="linear_trend", w=10) # A seasonality -seasonal52 = seasonality(period=52, discount=0.99, name='seasonal52', w=10) +seasonal52 = seasonality(period=52, discount=0.99, name="seasonal52", w=10) simple_dlm = dlm(time_series) + linear_trend + seasonal52 simple_dlm.fit() # Plot the fitted results -simple_dlm.turnOff('data points') +simple_dlm.turnOff("data points") simple_dlm.plot() # Plot each component (attribution) -simple_dlm.turnOff('predict plot') -simple_dlm.turnOff('filtered plot') -simple_dlm.plot('linear_trend') -simple_dlm.plot('seasonal52') +simple_dlm.turnOff("predict plot") +simple_dlm.turnOff("filtered plot") +simple_dlm.plot("linear_trend") +simple_dlm.plot("seasonal52") # Plot the prediction give the first 350 weeks and forcast the next 200 weeks. simple_dlm.plotPredictN(N=200, date=350) # Plot the prediction give the first 250 weeks and forcast the next 200 weeks. @@ -64,19 +67,20 @@ # Build a dynamic regression model from pydlm import dynamic -regressor10 = dynamic(features=features, discount=1.0, name='regressor10', w=10) + +regressor10 = dynamic(features=features, discount=1.0, name="regressor10", w=10) drm = dlm(time_series) + linear_trend + seasonal52 + regressor10 drm.fit() # Plot the fitted results -drm.turnOff('data points') +drm.turnOff("data points") drm.plot() # Plot each component (attribution) -drm.turnOff('predict plot') -drm.turnOff('filtered plot') -drm.plot('linear_trend') -drm.plot('seasonal52') -drm.plot('regressor10') +drm.turnOff("predict plot") +drm.turnOff("filtered plot") +drm.plot("linear_trend") +drm.plot("seasonal52") +drm.plot("regressor10") # Plot the prediction give the first 300 weeks and forcast the next 150 weeks. drm.plotPredictN(N=150, date=300) # Plot the prediction give the first 250 weeks and forcast the next 200 weeks. diff --git a/pydlm/__init__.py b/pydlm/__init__.py index 84282a72..fc4d750b 100644 --- a/pydlm/__init__.py +++ b/pydlm/__init__.py @@ -1,6 +1,14 @@ # This is the PyDLM package -__all__ = ['dlm', 'trend', 'seasonality', 'dynamic', 'autoReg', 'longSeason', 'modelTuner'] +__all__ = [ + "dlm", + "trend", + "seasonality", + "dynamic", + "autoReg", + "longSeason", + "modelTuner", +] from pydlm.dlm import dlm from pydlm.modeler.trends import trend diff --git a/pydlm/access/_dlmGet.py b/pydlm/access/_dlmGet.py index 8791cc75..fb63d389 100644 --- a/pydlm/access/_dlmGet.py +++ b/pydlm/access/_dlmGet.py @@ -6,12 +6,13 @@ =============================================================================== """ + from numpy import dot from pydlm.core._dlm import _dlm class _dlmGet(_dlm): - """ The class containing all get methods for dlm class. + """The class containing all get methods for dlm class. Methods: _getComponent: get the component if it is in dlm @@ -21,10 +22,9 @@ class _dlmGet(_dlm): _getComponentVar: get the variance of a given component """ - # function to get the corresponding latent state def _getLatentState(self, name, filterType, start, end): - """ Get the latent states of a given component. + """Get the latent states of a given component. Args: name: the name of the component. @@ -39,21 +39,20 @@ def _getLatentState(self, name, filterType, start, end): """ end += 1 indx = self.builder.componentIndex[name] - patten = lambda x: x if x is None else x[indx[0]:(indx[1] + 1), 0:1] + patten = lambda x: x if x is None else x[indx[0] : (indx[1] + 1), 0:1] - if filterType == 'forwardFilter': + if filterType == "forwardFilter": return list(map(patten, self.result.filteredState[start:end])) - elif filterType == 'backwardSmoother': + elif filterType == "backwardSmoother": return list(map(patten, self.result.smoothedState[start:end])) - elif filterType == 'predict': + elif filterType == "predict": return list(map(patten, self.result.predictedState[start:end])) else: - raise NameError('Incorrect filter type') - + raise NameError("Incorrect filter type") # function to get the corresponding latent covariance def _getLatentCov(self, name, filterType, start, end): - """ Get the latent covariance of a given component. + """Get the latent covariance of a given component. Args: name: the name of the component. @@ -68,22 +67,24 @@ def _getLatentCov(self, name, filterType, start, end): """ end += 1 indx = self.builder.componentIndex[name] - patten = lambda x: x if x is None \ - else x[indx[0]:(indx[1] + 1), indx[0]:(indx[1] + 1)] + patten = ( + lambda x: x + if x is None + else x[indx[0] : (indx[1] + 1), indx[0] : (indx[1] + 1)] + ) - if filterType == 'forwardFilter': + if filterType == "forwardFilter": return list(map(patten, self.result.filteredCov[start:end])) - elif filterType == 'backwardSmoother': + elif filterType == "backwardSmoother": return list(map(patten, self.result.smoothedCov[start:end])) - elif filterType == 'predict': + elif filterType == "predict": return list(map(patten, self.result.predictedCov[start:end])) else: - raise NameError('Incorrect filter type') - + raise NameError("Incorrect filter type") # function to get the component mean def _getComponentMean(self, name, filterType, start, end): - """ Get the mean of a given component. + """Get the mean of a given component. Args: name: the name of the component. @@ -98,23 +99,21 @@ def _getComponentMean(self, name, filterType, start, end): """ end += 1 comp = self._fetchComponent(name) - componentState = self._getLatentState(name=name, - filterType=filterType, - start=start, end=end) + componentState = self._getLatentState( + name=name, filterType=filterType, start=start, end=end + ) result = [] for k, i in enumerate(range(start, end)): if name in self.builder.dynamicComponents: comp.updateEvaluation(i) elif name in self.builder.automaticComponents: comp.updateEvaluation(i, self.padded_data) - result.append(dot(comp.evaluation, - componentState[k]).tolist()[0][0]) + result.append(dot(comp.evaluation, componentState[k]).tolist()[0][0]) return result - # function to get the component variance def _getComponentVar(self, name, filterType, start, end): - """ Get the variance of a given component. + """Get the variance of a given component. Args: name: the name of the component. @@ -129,16 +128,18 @@ def _getComponentVar(self, name, filterType, start, end): """ end += 1 comp = self._fetchComponent(name) - componentCov = self._getLatentCov(name=name, - filterType=filterType, - start=start, end=end) + componentCov = self._getLatentCov( + name=name, filterType=filterType, start=start, end=end + ) result = [] for k, i in enumerate(range(start, end)): if name in self.builder.dynamicComponents: comp.updateEvaluation(i) elif name in self.builder.automaticComponents: comp.updateEvaluation(i, self.padded_data) - result.append(dot( - dot(comp.evaluation, - componentCov[k]), comp.evaluation.T).tolist()[0][0]) + result.append( + dot(dot(comp.evaluation, componentCov[k]), comp.evaluation.T).tolist()[ + 0 + ][0] + ) return result diff --git a/pydlm/access/dlmAccessMod.py b/pydlm/access/dlmAccessMod.py index cb1ded4e..11c11d0c 100644 --- a/pydlm/access/dlmAccessMod.py +++ b/pydlm/access/dlmAccessMod.py @@ -4,14 +4,13 @@ class dlmAccessModule(_dlmGet): - """ A dlm module for all the access methods. This is an API layer for the + """A dlm module for all the access methods. This is an API layer for the `dlm` class. All methods defined in this class are public and can be called directly from `dlm` object. """ - def getAll(self): - """ get all the _result class which contains all results + """get all the _result class which contains all results Returns: The @result object containing all computed results. @@ -19,9 +18,8 @@ def getAll(self): """ return deepcopy(self.result) - - def getMean(self, filterType='forwardFilter', name='main'): - """ get mean for data or component. + def getMean(self, filterType="forwardFilter", name="main"): + """get mean for data or component. If the working dates are not (0, self.n - 1), then a warning will prompt stating the actual filtered dates. @@ -41,31 +39,27 @@ def getMean(self, filterType='forwardFilter', name='main'): """ # get the working date start, end = self._checkAndGetWorkingDates(filterType=filterType) - end += 1 # To get the result for the last date. + end += 1 # To get the result for the last date. # get the mean for the fitlered data - if name == 'main': + if name == "main": # get out of the matrix form - if filterType == 'forwardFilter': - return self._1DmatrixToArray( - self.result.filteredObs[start:end]) - elif filterType == 'backwardSmoother': - return self._1DmatrixToArray( - self.result.smoothedObs[start:end]) - elif filterType == 'predict': - return self._1DmatrixToArray( - self.result.predictedObs[start:end]) + if filterType == "forwardFilter": + return self._1DmatrixToArray(self.result.filteredObs[start:end]) + elif filterType == "backwardSmoother": + return self._1DmatrixToArray(self.result.smoothedObs[start:end]) + elif filterType == "predict": + return self._1DmatrixToArray(self.result.predictedObs[start:end]) else: - raise NameError('Incorrect filter type.') + raise NameError("Incorrect filter type.") # get the mean for the component self._checkComponent(name) - return self._getComponentMean(name=name, - filterType=filterType, - start=start, end=(end - 1)) + return self._getComponentMean( + name=name, filterType=filterType, start=start, end=(end - 1) + ) - - def getVar(self, filterType='forwardFilter', name='main'): - """ get the variance for data or component. + def getVar(self, filterType="forwardFilter", name="main"): + """get the variance for data or component. If the filtered dates are not (0, self.n - 1), then a warning will prompt stating @@ -88,28 +82,25 @@ def getVar(self, filterType='forwardFilter', name='main'): start, end = self._checkAndGetWorkingDates(filterType=filterType) end += 1 # get the variance for the time series data - if name == 'main': + if name == "main": # get out of the matrix form - if filterType == 'forwardFilter': - return self._1DmatrixToArray( - self.result.filteredObsVar[start:end]) - elif filterType == 'backwardSmoother': - return self._1DmatrixToArray( - self.result.smoothedObsVar[start:end]) - elif filterType == 'predict': - return self._1DmatrixToArray( - self.result.predictedObsVar[start:end]) + if filterType == "forwardFilter": + return self._1DmatrixToArray(self.result.filteredObsVar[start:end]) + elif filterType == "backwardSmoother": + return self._1DmatrixToArray(self.result.smoothedObsVar[start:end]) + elif filterType == "predict": + return self._1DmatrixToArray(self.result.predictedObsVar[start:end]) else: - raise NameError('Incorrect filter type.') + raise NameError("Incorrect filter type.") # get the variance for the component self._checkComponent(name) - return self._getComponentVar(name=name, filterType=filterType, - start=start, end=(end - 1)) - + return self._getComponentVar( + name=name, filterType=filterType, start=start, end=(end - 1) + ) - def getResidual(self, filterType='forwardFilter'): - """ get the residuals for data after filtering or smoothing. + def getResidual(self, filterType="forwardFilter"): + """get the residuals for data after filtering or smoothing. If the working dates are not (0, self.n - 1), then a warning will prompt stating the actual filtered dates. @@ -125,27 +116,26 @@ def getResidual(self, filterType='forwardFilter'): """ # get the working date start, end = self._checkAndGetWorkingDates(filterType=filterType) - end += 1 # To get the result for the last date. + end += 1 # To get the result for the last date. # get the mean for the fitlered data # get out of the matrix form - if filterType == 'forwardFilter': + if filterType == "forwardFilter": return self._1DmatrixToArray( - [self.data[i] - self.result.filteredObs[i] - for i in range(start, end)]) - elif filterType == 'backwardSmoother': + [self.data[i] - self.result.filteredObs[i] for i in range(start, end)] + ) + elif filterType == "backwardSmoother": return self._1DmatrixToArray( - [self.data[i] - self.result.smoothedObs[i] - for i in range(start, end)]) - elif filterType == 'predict': + [self.data[i] - self.result.smoothedObs[i] for i in range(start, end)] + ) + elif filterType == "predict": return self._1DmatrixToArray( - [self.data[i] - self.result.predictedObs[i] - for i in range(start, end)]) + [self.data[i] - self.result.predictedObs[i] for i in range(start, end)] + ) else: - raise NameError('Incorrect filter type.') + raise NameError("Incorrect filter type.") - - def getInterval(self, p=0.95, filterType='forwardFilter', name='main'): - """ get the confidence interval for data or component. + def getInterval(self, p=0.95, filterType="forwardFilter", name="main"): + """get the confidence interval for data or component. If the filtered dates are not (0, self.n - 1), then a warning will prompt stating the actual @@ -170,43 +160,36 @@ def getInterval(self, p=0.95, filterType='forwardFilter', name='main'): start, end = self._checkAndGetWorkingDates(filterType=filterType) end += 1 # get the mean and the variance for the time series data - if name == 'main': + if name == "main": # get out of the matrix form - if filterType == 'forwardFilter': - compMean = self._1DmatrixToArray( - self.result.filteredObs[start:end]) - compVar = self._1DmatrixToArray( - self.result.filteredObsVar[start:end]) - elif filterType == 'backwardSmoother': - compMean = self._1DmatrixToArray( - self.result.smoothedObs[start:end]) - compVar = self._1DmatrixToArray( - self.result.smoothedObsVar[start:end]) - elif filterType == 'predict': - compMean = self._1DmatrixToArray( - self.result.predictedObs[start:end]) - compVar = self._1DmatrixToArray( - self.result.predictedObsVar[start:end]) + if filterType == "forwardFilter": + compMean = self._1DmatrixToArray(self.result.filteredObs[start:end]) + compVar = self._1DmatrixToArray(self.result.filteredObsVar[start:end]) + elif filterType == "backwardSmoother": + compMean = self._1DmatrixToArray(self.result.smoothedObs[start:end]) + compVar = self._1DmatrixToArray(self.result.smoothedObsVar[start:end]) + elif filterType == "predict": + compMean = self._1DmatrixToArray(self.result.predictedObs[start:end]) + compVar = self._1DmatrixToArray(self.result.predictedObsVar[start:end]) else: - raise NameError('Incorrect filter type.') + raise NameError("Incorrect filter type.") # get the mean and variance for the component else: self._checkComponent(name) - compMean = self._getComponentMean(name=name, - filterType=filterType, - start=start, end=(end - 1)) - compVar = self._getComponentVar(name=name, - filterType=filterType, - start=start, end=(end - 1)) + compMean = self._getComponentMean( + name=name, filterType=filterType, start=start, end=(end - 1) + ) + compVar = self._getComponentVar( + name=name, filterType=filterType, start=start, end=(end - 1) + ) # get the upper and lower bound upper, lower = getInterval(compMean, compVar, p) return (upper, lower) - - def getLatentState(self, filterType='forwardFilter', name='all'): - """ get the latent states for different components and filters. + def getLatentState(self, filterType="forwardFilter", name="all"): + """get the latent states for different components and filters. If the filtered dates are not (0, self.n - 1), then a warning will prompt stating the actual filtered dates. @@ -230,24 +213,29 @@ def getLatentState(self, filterType='forwardFilter', name='all'): end += 1 # to return the full latent states to_array = lambda x: None if x is None else x.flatten().tolist() - if name == 'all': - if filterType == 'forwardFilter': + if name == "all": + if filterType == "forwardFilter": return list(map(to_array, self.result.filteredState[start:end])) - elif filterType == 'backwardSmoother': + elif filterType == "backwardSmoother": return list(map(to_array, self.result.smoothedState[start:end])) - elif filterType == 'predict': + elif filterType == "predict": return list(map(to_array, self.result.predictedState[start:end])) else: - raise NameError('Incorrect filter type.') + raise NameError("Incorrect filter type.") # to return the latent state for a given component self._checkComponent(name) - return list(map(to_array, self._getLatentState(name=name, filterType=filterType, - start=start, end=(end - 1)))) - - - def getLatentCov(self, filterType='forwardFilter', name='all'): - """ get the error covariance for different components and + return list( + map( + to_array, + self._getLatentState( + name=name, filterType=filterType, start=start, end=(end - 1) + ), + ) + ) + + def getLatentCov(self, filterType="forwardFilter", name="all"): + """get the error covariance for different components and filters. If the filtered dates are not (0, self.n - 1), @@ -271,17 +259,18 @@ def getLatentCov(self, filterType='forwardFilter', name='all'): start, end = self._checkAndGetWorkingDates(filterType=filterType) end += 1 # to return the full latent covariance - if name == 'all': - if filterType == 'forwardFilter': + if name == "all": + if filterType == "forwardFilter": return self.result.filteredCov[start:end] - elif filterType == 'backwardSmoother': + elif filterType == "backwardSmoother": return self.result.smoothedCov[start:end] - elif filterType == 'predict': + elif filterType == "predict": return self.result.predictedCov[start:end] else: - raise NameError('Incorrect filter type.') + raise NameError("Incorrect filter type.") # to return the latent covariance for a given component self._checkComponent(name) - return self._getLatentCov(name=name, filterType=filterType, - start=start, end=(end - 1)) + return self._getLatentCov( + name=name, filterType=filterType, start=start, end=(end - 1) + ) diff --git a/pydlm/base/__init__.py b/pydlm/base/__init__.py index 993f3c7c..0cee8188 100644 --- a/pydlm/base/__init__.py +++ b/pydlm/base/__init__.py @@ -1,7 +1,7 @@ # the base kalmanFilter for filtering -#__all__ = ['kalmanFilter', 'baseModel', 'tools'] +# __all__ = ['kalmanFilter', 'baseModel', 'tools'] -#import pydlm.base.kalmanFilter -#from .tools import tools -#from .baseModel import baseModel +# import pydlm.base.kalmanFilter +# from .tools import tools +# from .baseModel import baseModel diff --git a/pydlm/base/baseModel.py b/pydlm/base/baseModel.py index ef71d6f4..0f0949fc 100644 --- a/pydlm/base/baseModel.py +++ b/pydlm/base/baseModel.py @@ -9,6 +9,7 @@ It stores all the necessary components for kalmanFilter and save the results """ + # dependencies import numpy as np import pydlm.base.tools as tl @@ -16,7 +17,7 @@ # define the basic structure for a dlm model class baseModel: - """ The baseModel class that provides the basic model structure for dlm. + """The baseModel class that provides the basic model structure for dlm. Attributes: transition: the transition matrix G @@ -34,10 +35,17 @@ class baseModel: validation: validate the matrix dimensions are consistent. """ - # define the components of a baseModel - def __init__(self, transition = None, evaluation = None, noiseVar = None, \ - sysVar = None, innovation = None, state = None, df = None): + def __init__( + self, + transition=None, + evaluation=None, + noiseVar=None, + sysVar=None, + innovation=None, + state=None, + df=None, + ): self.transition = transition self.evaluation = evaluation self.noiseVar = noiseVar @@ -51,23 +59,19 @@ def __init__(self, transition = None, evaluation = None, noiseVar = None, \ # a hidden data field used only for model prediction self.prediction = __model__() - # initialize the observation mean and variance def initializeObservation(self): - """ Initialize the value of obs and obsVar - - """ + """Initialize the value of obs and obsVar""" self.validation() self.obs = np.dot(self.evaluation, self.state) - self.obsVar = np.dot(np.dot(self.evaluation, self.sysVar), self.evaluation.T) \ - + self.noiseVar - + self.obsVar = ( + np.dot(np.dot(self.evaluation, self.sysVar), self.evaluation.T) + + self.noiseVar + ) # checking if the dimension matches with each other def validation(self): - """ Validate the model components are consistent - - """ + """Validate the model components are consistent""" # check symmetric tl.checker.checkSymmetry(self.transition) tl.checker.checkSymmetry(self.sysVar) @@ -81,12 +85,11 @@ def validation(self): tl.checker.checkVectorDimension(self.evaluation, self.transition) tl.checker.checkVectorDimension(self.state, self.transition) - + # define an inner class to store intermediate results class __model__: - # to store result for prediction - def __init__(self, step = 0, state = None, obs = None, sysVar = None, obsVar = None): + def __init__(self, step=0, state=None, obs=None, sysVar=None, obsVar=None): self.step = 0 self.state = state self.obs = obs diff --git a/pydlm/base/kalmanFilter.py b/pydlm/base/kalmanFilter.py index 43805702..4f207e06 100644 --- a/pydlm/base/kalmanFilter.py +++ b/pydlm/base/kalmanFilter.py @@ -12,6 +12,7 @@ provided) """ + # This code take care of the Kalman filter import numpy as np import pydlm.base.tools as tl @@ -21,7 +22,7 @@ class kalmanFilter: - """ The kalmanFilter class the provide the basic functionalities + """The kalmanFilter class the provide the basic functionalities Attributes: discount: the discounting factor determining how much information to carry on @@ -38,11 +39,8 @@ class kalmanFilter: updateDiscount: for updating the discount factors """ - - def __init__(self, discount=[0.99], \ - updateInnovation='whole', - index=None): - """ Initializing the kalmanFilter class + def __init__(self, discount=[0.99], updateInnovation="whole", index=None): + """Initializing the kalmanFilter class Args: discount: the discounting factor, could be a vector @@ -55,9 +53,8 @@ def __init__(self, discount=[0.99], \ self.updateInnovation = updateInnovation self.index = index - - def predict(self, model, dealWithMissingEvaluation = False): - """ Predict the next states of the model by one step + def predict(self, model, dealWithMissingEvaluation=False): + """Predict the next states of the model by one step Args: model: the @baseModel class provided all necessary information @@ -74,45 +71,52 @@ def predict(self, model, dealWithMissingEvaluation = False): # if the step number == 0, we use result from the model state if model.prediction.step == 0: - model.prediction.state = np.dot(model.transition, model.state) model.prediction.obs = np.dot(model.evaluation, model.prediction.state) - model.prediction.sysVar = np.dot(np.dot(model.transition, model.sysVar), - model.transition.T) + model.prediction.sysVar = np.dot( + np.dot(model.transition, model.sysVar), model.transition.T + ) # update the innovation - if self.updateInnovation == 'whole': + if self.updateInnovation == "whole": self.__updateInnovation__(model) - elif self.updateInnovation == 'component': + elif self.updateInnovation == "component": self.__updateInnovation2__(model) # add the innovation to the system variance model.prediction.sysVar += model.innovation - model.prediction.obsVar = np.dot(np.dot(model.evaluation, - model.prediction.sysVar), - model.evaluation.T) + model.noiseVar + model.prediction.obsVar = ( + np.dot( + np.dot(model.evaluation, model.prediction.sysVar), + model.evaluation.T, + ) + + model.noiseVar + ) model.prediction.step = 1 # otherwise, we use previous result to predict next time stamp else: model.prediction.state = np.dot(model.transition, model.prediction.state) model.prediction.obs = np.dot(model.evaluation, model.prediction.state) - model.prediction.sysVar = np.dot(np.dot(model.transition, - model.prediction.sysVar), - model.transition.T) - model.prediction.obsVar = np.dot(np.dot(model.evaluation, - model.prediction.sysVar), - model.evaluation.T) + model.noiseVar + model.prediction.sysVar = np.dot( + np.dot(model.transition, model.prediction.sysVar), model.transition.T + ) + model.prediction.obsVar = ( + np.dot( + np.dot(model.evaluation, model.prediction.sysVar), + model.evaluation.T, + ) + + model.noiseVar + ) model.prediction.step += 1 # recover the evaluation and the transition matrix if dealWithMissingEvaluation: self._recoverTransitionAndEvaluation(model, loc) - - def forwardFilter(self, model, y, dealWithMissingEvaluation = False): - """ The forwardFilter used to run one step filtering given new data + def forwardFilter(self, model, y, dealWithMissingEvaluation=False): + """The forwardFilter used to run one step filtering given new data Args: model: the @baseModel provided the basic information @@ -131,7 +135,6 @@ def forwardFilter(self, model, y, dealWithMissingEvaluation = False): # when y is not a missing data if y is not None: - # first obtain the predicted status # we make the prediction step equal to 0 to ensure the prediction # is based on the model state and innovation is updated correctlly @@ -140,25 +143,34 @@ def forwardFilter(self, model, y, dealWithMissingEvaluation = False): # the prediction error and the correction matrix err = y - model.prediction.obs - correction = (np.dot(model.prediction.sysVar, model.evaluation.T) - / model.prediction.obsVar) + correction = ( + np.dot(model.prediction.sysVar, model.evaluation.T) + / model.prediction.obsVar + ) # update new states model.df += 1 - lastNoiseVar = model.noiseVar # for updating model.sysVar - model.noiseVar = (model.noiseVar * - (1.0 - 1.0 / model.df + - err * err / model.df / model.prediction.obsVar)) + lastNoiseVar = model.noiseVar # for updating model.sysVar + model.noiseVar = model.noiseVar * ( + 1.0 - 1.0 / model.df + err * err / model.df / model.prediction.obsVar + ) model.state = model.prediction.state + correction * err - model.sysVar = (model.noiseVar[0, 0] / lastNoiseVar[0, 0] * - (model.prediction.sysVar - np.dot(correction, correction.T) * - model.prediction.obsVar[0, 0])) + model.sysVar = ( + model.noiseVar[0, 0] + / lastNoiseVar[0, 0] + * ( + model.prediction.sysVar + - np.dot(correction, correction.T) * model.prediction.obsVar[0, 0] + ) + ) model.obs = np.dot(model.evaluation, model.state) - model.obsVar = np.dot(np.dot(model.evaluation, model.sysVar), - model.evaluation.T) + model.noiseVar + model.obsVar = ( + np.dot(np.dot(model.evaluation, model.sysVar), model.evaluation.T) + + model.noiseVar + ) # update the innovation using discount # model.innovation = model.sysVar * (1 / self.discount - 1) @@ -183,7 +195,6 @@ def forwardFilter(self, model, y, dealWithMissingEvaluation = False): if dealWithMissingEvaluation: self._recoverTransitionAndEvaluation(model, loc) - # The backward smoother for a given unsmoothed states at time t # what model should store: # model.state: the last smoothed states (t + 1) @@ -195,8 +206,7 @@ def forwardFilter(self, model, y, dealWithMissingEvaluation = False): # rawState: the unsmoothed state at time t # rawSysVar: the unsmoothed system variance at time t def backwardSmoother(self, model, rawState, rawSysVar): - - """ The backwardSmoother for one step backward smoothing + """The backwardSmoother for one step backward smoothing Args: model: the @baseModel used for backward smoothing, the model shall store @@ -227,21 +237,24 @@ def backwardSmoother(self, model, rawState, rawSysVar): ################################################ backward = np.dot(np.dot(rawSysVar, model.transition.T), predSysVarInv) - model.state = rawState + np.dot(backward, (model.state - model.prediction.state)) - model.sysVar = (rawSysVar + - np.dot(np.dot(backward, - (model.sysVar - model.prediction.sysVar)), backward.T)) + model.state = rawState + np.dot( + backward, (model.state - model.prediction.state) + ) + model.sysVar = rawSysVar + np.dot( + np.dot(backward, (model.sysVar - model.prediction.sysVar)), backward.T + ) model.obs = np.dot(model.evaluation, model.state) - model.obsVar = np.dot(np.dot(model.evaluation, model.sysVar), - model.evaluation.T) + model.noiseVar + model.obsVar = ( + np.dot(np.dot(model.evaluation, model.sysVar), model.evaluation.T) + + model.noiseVar + ) # recover the evaluation and the transition matrix - #if dealWithMissingEvaluation: + # if dealWithMissingEvaluation: # self._recoverTransitionAndEvaluation(model, loc) - def backwardSampler(self, model, rawState, rawSysVar): - """ The backwardSampler for one step backward sampling + """The backwardSampler for one step backward sampling Args: model: the @baseModel used for backward sampling, the model shall store @@ -267,22 +280,23 @@ def backwardSampler(self, model, rawState, rawSysVar): ################################################ backward = np.dot(np.dot(rawSysVar, model.transition.T), predSysVarInv) - model.state = rawState + np.dot(backward, (model.state - model.prediction.state)) - model.sysVar = (rawSysVar + - np.dot(np.dot(backward, - (model.sysVar - model.prediction.sysVar)), backward.T)) - model.state = np.random.multivariate_normal(model.state.A1, - model.sysVar).T + model.state = rawState + np.dot( + backward, (model.state - model.prediction.state) + ) + model.sysVar = rawSysVar + np.dot( + np.dot(backward, (model.sysVar - model.prediction.sysVar)), backward.T + ) + model.state = np.random.multivariate_normal(model.state.A1, model.sysVar).T model.obs = np.dot(model.evaluation, model.state) - model.obsVar = np.dot(np.dot(model.evaluation, model.sysVar), - model.evaluation.T) + model.noiseVar - model.obs = np.random.multivariate_normal(model.obs.A1, - model.obsVar).T - + model.obsVar = ( + np.dot(np.dot(model.evaluation, model.sysVar), model.evaluation.T) + + model.noiseVar + ) + model.obs = np.random.multivariate_normal(model.obs.A1, model.obsVar).T # for updating the discounting factor def updateDiscount(self, newDiscount): - """ For updating the discounting factor + """For updating the discounting factor Args: newDiscount: the new discount factor @@ -291,54 +305,48 @@ def updateDiscount(self, newDiscount): self.__checkDiscount__(newDiscount) self.discount = np.diag(1 / np.sqrt(newDiscount)) - def __checkDiscount__(self, discount): - """ Check whether the discount fact is within (0, 1) - - """ + """Check whether the discount fact is within (0, 1)""" for i in range(len(discount)): if discount[i] < 0 or discount[i] > 1: - raise tl.matrixErrors('discount factor must be between 0 and 1') - + raise tl.matrixErrors("discount factor must be between 0 and 1") # update the innovation def __updateInnovation__(self, model): - """ update the innovation matrix of the model - - """ - - model.innovation = np.dot(np.dot(self.discount, model.prediction.sysVar), - self.discount) - model.prediction.sysVar + """update the innovation matrix of the model""" + model.innovation = ( + np.dot(np.dot(self.discount, model.prediction.sysVar), self.discount) + - model.prediction.sysVar + ) # update the innovation def __updateInnovation2__(self, model): - """ update the innovation matrix of the model, but only for component + """update the innovation matrix of the model, but only for component indepdently. (i.e., only add innovation to block diagonals, not on off block diagonals) """ - innovation = np.dot(np.dot(self.discount, model.prediction.sysVar), - self.discount) - model.prediction.sysVar + innovation = ( + np.dot(np.dot(self.discount, model.prediction.sysVar), self.discount) + - model.prediction.sysVar + ) model.innovation = np.zeros(innovation.shape) for name in self.index: indx = self.index[name] - model.innovation[indx[0]: (indx[1] + 1), indx[0]: (indx[1] + 1)] = ( - innovation[indx[0]: (indx[1] + 1), indx[0]: (indx[1] + 1)]) - + model.innovation[indx[0] : (indx[1] + 1), indx[0] : (indx[1] + 1)] = ( + innovation[indx[0] : (indx[1] + 1), indx[0] : (indx[1] + 1)] + ) # a generalized inverse of matrix A def _gInverse(self, A): - """ A generalized inverse of matrix A - - """ + """A generalized inverse of matrix A""" return np.linalg.pinv(A) - def _modifyTransitionAccordingToMissingValue(self, model): - """ When evaluation contains None value, we modify the corresponding entries + """When evaluation contains None value, we modify the corresponding entries in the transition to deal with the missing value """ @@ -350,13 +358,11 @@ def _modifyTransitionAccordingToMissingValue(self, model): model.evaluation[0, i] = 0.0 return loc - def _recoverTransitionAndEvaluation(self, model, loc): - """ We recover the transition and evaluation use the results from + """We recover the transition and evaluation use the results from _modifyTransitionAccordingToMissingValue """ for i in loc: model.evaluation[0, i] = None model.transition[i, i] = 1.0 - diff --git a/pydlm/base/tools.py b/pydlm/base/tools.py index d9a42a33..dc82e7e7 100644 --- a/pydlm/base/tools.py +++ b/pydlm/base/tools.py @@ -3,7 +3,6 @@ # define the error class for exceptions class matrixErrors(Exception): - def __init__(self, value): self.value = value @@ -13,33 +12,33 @@ def __str__(self): # The class to check matrixErrors class checker: - - # checking if two matrix has the same dimension @staticmethod def checkMatrixDimension(A, B): if A.shape != B.shape: - raise matrixErrors('The dimensions do not match!') - + raise matrixErrors("The dimensions do not match!") # checking if a vector has the dimension as matrix @staticmethod def checkVectorDimension(v, A): vDim = v.shape ADim = A.shape - if vDim[0] != ADim[0] and vDim[0] != ADim[1] and \ - vDim[1] != ADim[0] and vDim[1] != ADim[1]: - raise matrixErrors('The dimensions do not match!') - + if ( + vDim[0] != ADim[0] + and vDim[0] != ADim[1] + and vDim[1] != ADim[0] + and vDim[1] != ADim[1] + ): + raise matrixErrors("The dimensions do not match!") # checking if a matrix is symmetric @staticmethod def checkSymmetry(A): ADim = A.shape if ADim[0] != ADim[1]: - raise matrixErrors('The matrix is not symmetric!') + raise matrixErrors("The matrix is not symmetric!") + - # a completely unshared copy of lists def duplicateList(aList): if isinstance(aList, list): @@ -49,26 +48,25 @@ def duplicateList(aList): # inverse normal cdf function def rational_approximation(t): - # Abramowitz and Stegun formula 26.2.23. # The absolute value of the error should be less than 4.5 e-4. c = [2.515517, 0.802853, 0.010328] d = [1.432788, 0.189269, 0.001308] - numerator = (c[2]*t + c[1])*t + c[0] - denominator = ((d[2]*t + d[1])*t + d[0])*t + 1.0 + numerator = (c[2] * t + c[1]) * t + c[0] + denominator = ((d[2] * t + d[1]) * t + d[0]) * t + 1.0 return t - numerator / denominator - - + + def normal_CDF_inverse(p): assert p > 0.0 and p < 1 - + # See article above for explanation of this section. if p < 0.5: # F^-1(p) = - G^-1(p) - return -rational_approximation( math.sqrt(-2.0*math.log(p)) ) + return -rational_approximation(math.sqrt(-2.0 * math.log(p))) else: # F^-1(p) = G^-1(1-p) - return rational_approximation( math.sqrt(-2.0*math.log(1.0-p)) ) + return rational_approximation(math.sqrt(-2.0 * math.log(1.0 - p))) def getInterval(means, var, p): diff --git a/pydlm/core/__init__.py b/pydlm/core/__init__.py index aa8b129b..119921b6 100644 --- a/pydlm/core/__init__.py +++ b/pydlm/core/__init__.py @@ -1,3 +1,3 @@ -#__all__ = ['_dlm'] +# __all__ = ['_dlm'] -#import pydlm.func._dlm as _dlm +# import pydlm.func._dlm as _dlm diff --git a/pydlm/core/_dlm.py b/pydlm/core/_dlm.py index b50fa9b5..a2ff07ec 100644 --- a/pydlm/core/_dlm.py +++ b/pydlm/core/_dlm.py @@ -9,6 +9,7 @@ It provides the basic modeling, filtering, forecasting and smoothing of a dlm. """ + from pydlm.base.kalmanFilter import kalmanFilter from pydlm.modeler.builder import builder @@ -21,7 +22,7 @@ class _dlm(object): - """ _dlm includes all basic functions for building dlm model. + """_dlm includes all basic functions for building dlm model. Attributes: data: the observed time series data @@ -58,11 +59,9 @@ class _dlm(object): _checkAndGetWorkingDates: get the correct filtering dates """ - # define the basic members # initialize the result def __init__(self, data, **options): - self.data = list(data) # padded_data is used by auto regressor. It is the raw data with missing value # replaced by forward filter results. (Missing value include the out of scope @@ -77,63 +76,65 @@ def __init__(self, data, **options): self.initialized = False self.time = None - # an inner class to store all options class _defaultOptions(object): - """ All plotting and fitting options + """All plotting and fitting options""" - """ def __init__(self, **kwargs): - self.noise = kwargs.get('noise', 1.0) + self.noise = kwargs.get("noise", 1.0) # Stable mode is basically doing rolling window refitting. Use with caution. - self.stable = kwargs.get('stable', False) - self.innovationType = kwargs.get('component', 'component') - - self.plotOriginalData = kwargs.get('plotOriginalData', True) - self.plotFilteredData = kwargs.get('plotFilteredData', True) - self.plotSmoothedData = kwargs.get('plotSmoothedData', True) - self.plotPredictedData = kwargs.get('plotPredictedData', True) - self.showDataPoint = kwargs.get('showDataPoint', False) - self.showFittedPoint = kwargs.get('showFittedPoint', False) - self.showConfidenceInterval = kwargs.get('showConfidenceInterval', True) - self.dataColor = kwargs.get('dataColor', 'black') - self.filteredColor = kwargs.get('filteredColor', 'blue') - self.predictedColor = kwargs.get('predictedColor', 'green') - self.smoothedColor = kwargs.get('smoothedColor', 'red') - self.separatePlot = kwargs.get('seperatePlot', True) - self.confidence = kwargs.get('confidence', 0.95) - self.intervalType = kwargs.get('intervalType', 'ribbon') - self.useAutoNoise = kwargs.get('useAutoNoise', False) - self.logger = kwargs.get('logger', self._getDefaultLogger()) + self.stable = kwargs.get("stable", False) + self.innovationType = kwargs.get("component", "component") + + self.plotOriginalData = kwargs.get("plotOriginalData", True) + self.plotFilteredData = kwargs.get("plotFilteredData", True) + self.plotSmoothedData = kwargs.get("plotSmoothedData", True) + self.plotPredictedData = kwargs.get("plotPredictedData", True) + self.showDataPoint = kwargs.get("showDataPoint", False) + self.showFittedPoint = kwargs.get("showFittedPoint", False) + self.showConfidenceInterval = kwargs.get("showConfidenceInterval", True) + self.dataColor = kwargs.get("dataColor", "black") + self.filteredColor = kwargs.get("filteredColor", "blue") + self.predictedColor = kwargs.get("predictedColor", "green") + self.smoothedColor = kwargs.get("smoothedColor", "red") + self.separatePlot = kwargs.get("seperatePlot", True) + self.confidence = kwargs.get("confidence", 0.95) + self.intervalType = kwargs.get("intervalType", "ribbon") + self.useAutoNoise = kwargs.get("useAutoNoise", False) + self.logger = kwargs.get("logger", self._getDefaultLogger()) def _getDefaultLogger(self): - """ Get default logger """ - + """Get default logger""" + logging.basicConfig() - logger = logging.getLogger('pydlm') + logger = logging.getLogger("pydlm") logger.setLevel(logging.INFO) return logger - - - # an inner class to store all results class _result(object): - """ Class to store the results + """Class to store the results""" - """ # class level (static) variables to record all names - records = ['filteredObs', 'predictedObs', 'smoothedObs', - 'filteredObsVar', - 'predictedObsVar', 'smoothedObsVar', 'noiseVar', - 'df', - 'filteredState', 'predictedState', 'smoothedState', - 'filteredCov', 'predictedCov', 'smoothedCov'] - + records = [ + "filteredObs", + "predictedObs", + "smoothedObs", + "filteredObsVar", + "predictedObsVar", + "smoothedObsVar", + "noiseVar", + "df", + "filteredState", + "predictedState", + "smoothedState", + "filteredCov", + "predictedCov", + "smoothedCov", + ] # quantites to record the result def __init__(self, n): - # initialize all records to be [None] * n for variable in self.records: setattr(self, variable, [None] * n) @@ -149,44 +150,43 @@ def __init__(self, n): self.predictStatus = None # One day ahead prediction squared error - # extend the current record by n blocks def _appendResult(self, n): for variable in self.records: getattr(self, variable).extend([None] * n) - # pop out a specific date def _popout(self, date): for variable in self.records: getattr(self, variable).pop(date) - # initialize the builder def _initialize(self): - """ Initialize the model: initialize builder and filter. - - """ + """Initialize the model: initialize builder and filter.""" self._autoNoise() self.builder.initialize(noise=self.options.noise, data=self.padded_data) - self.Filter = kalmanFilter(discount=self.builder.discount, - updateInnovation=self.options.innovationType, - index=self.builder.componentIndex) + self.Filter = kalmanFilter( + discount=self.builder.discount, + updateInnovation=self.options.innovationType, + index=self.builder.componentIndex, + ) self.result = self._result(self.n) self.initialized = True - def _initializeFromBuilder(self, exported_builder): - self.builder.initializeFromBuilder(data=self.data, exported_builder=exported_builder) - self.Filter = kalmanFilter(discount=self.builder.discount, - updateInnovation=self.options.innovationType, - index=self.builder.componentIndex) + self.builder.initializeFromBuilder( + data=self.data, exported_builder=exported_builder + ) + self.Filter = kalmanFilter( + discount=self.builder.discount, + updateInnovation=self.options.innovationType, + index=self.builder.componentIndex, + ) self.result = self._result(self.n) self.initialized = True - def _autoNoise(self): - """ Auto initialize the noise parameter if options.useAutoNoise + """Auto initialize the noise parameter if options.useAutoNoise is true. """ @@ -194,21 +194,17 @@ def _autoNoise(self): trimmed_data = [x for x in self.data if x is not None] self.options.noise = min(var(trimmed_data), 1) - - # use the forward filter to filter the data + # use the forward filter to filter the data # start: the place where the filter started # end: the place where the filter ended # save: the index for dates where the filtered results should be saved, # could be 'all' or 'end' # isForget: indicate where the filter should use the previous state as # prior or just use the prior information from builder - def _forwardFilter(self, - start=0, - end=None, - save='all', - ForgetPrevious=False, - renew=False): - """ Running forwardFilter for the data for a given start and end date + def _forwardFilter( + self, start=0, end=None, save="all", ForgetPrevious=False, renew=False + ): + """Running forwardFilter for the data for a given start and end date Args: start: the start date @@ -244,48 +240,53 @@ def _forwardFilter(self, # otherwise we use the result on the previous day as the prior else: if start > self.result.filteredSteps[1] + 1: - raise ValueError('The data before start date has' - ' yet to be filtered! Otherwise set' - ' ForgetPrevious to be True. Check the' - ' in object.') + raise ValueError( + "The data before start date has" + " yet to be filtered! Otherwise set" + " ForgetPrevious to be True. Check the" + " in object." + ) self._setModelStatus(date=start - 1) # we run the forward filter sequentially lastRenewPoint = start # record the last renew point for step in range(start, end + 1): - # first check whether we need to update evaluation or not - if len(self.builder.dynamicComponents) > 0 or \ - len(self.builder.automaticComponents) > 0: + if ( + len(self.builder.dynamicComponents) > 0 + or len(self.builder.automaticComponents) > 0 + ): self.builder.updateEvaluation(step, self.padded_data) # check if rewnew is needed - if renew and step - lastRenewPoint > self.builder.renewTerm \ - and self.builder.renewTerm > 0.0: + if ( + renew + and step - lastRenewPoint > self.builder.renewTerm + and self.builder.renewTerm > 0.0 + ): # we renew the state of the day self._resetModelStatus() - for innerStep in range(step - int(self.builder.renewTerm), - step): - self.Filter.forwardFilter(self.builder.model, - self.data[innerStep]) + for innerStep in range(step - int(self.builder.renewTerm), step): + self.Filter.forwardFilter(self.builder.model, self.data[innerStep]) lastRenewPoint = step # then we use the updated model to filter the state self.Filter.forwardFilter(self.builder.model, self.data[step]) # extract the result and record - if save == 'all' or save == step: - self._copy(model=self.builder.model, - result=self.result, - step=step, - filterType='forwardFilter') - + if save == "all" or save == step: + self._copy( + model=self.builder.model, + result=self.result, + step=step, + filterType="forwardFilter", + ) # use the backward smooth to smooth the state # start: the last date of the backward filtering chain # days: number of days to go back from start def _backwardSmoother(self, start=None, days=None, ignoreFuture=False): - """ Backward smooth over filtered results for a specific start + """Backward smooth over filtered results for a specific start and number of days Args: @@ -308,16 +309,17 @@ def _backwardSmoother(self, start=None, days=None, ignoreFuture=False): # the forwardFilter has to be run before the smoother if self.result.filteredSteps[1] < start: - raise valueError('The last day has to be filtered before smoothing!' - 'check the in object.') + raise valueError( + "The last day has to be filtered before smoothing!" + "check the in object." + ) # and we record the most recent day which does not need to be smooth if start == self.n - 1 or ignoreFuture is True: self.result.smoothedState[start] = self.result.filteredState[start] self.result.smoothedObs[start] = self.result.filteredObs[start] self.result.smoothedCov[start] = self.result.filteredCov[start] - self.result.smoothedObsVar[start] \ - = self.result.filteredObsVar[start] + self.result.smoothedObsVar[start] = self.result.filteredObsVar[start] self.builder.model.noiseVar = self.result.noiseVar[start] start -= 1 else: @@ -336,31 +338,33 @@ def _backwardSmoother(self, start=None, days=None, ignoreFuture=False): dates.reverse() for day in dates: # we first update the model to be correct status before smooth - self.builder.model.prediction.state \ - = self.result.predictedState[day + 1] - self.builder.model.prediction.sysVar \ - = self.result.predictedCov[day + 1] + self.builder.model.prediction.state = self.result.predictedState[day + 1] + self.builder.model.prediction.sysVar = self.result.predictedCov[day + 1] - if len(self.builder.dynamicComponents) > 0 or \ - len(self.builder.automaticComponents) > 0: + if ( + len(self.builder.dynamicComponents) > 0 + or len(self.builder.automaticComponents) > 0 + ): self.builder.updateEvaluation(day, self.padded_data) # then we use the backward filter to filter the result self.Filter.backwardSmoother( model=self.builder.model, rawState=self.result.filteredState[day], - rawSysVar=self.result.filteredCov[day]) + rawSysVar=self.result.filteredCov[day], + ) # extract the result - self._copy(model=self.builder.model, - result=self.result, - step=day, - filterType='backwardSmoother') - + self._copy( + model=self.builder.model, + result=self.result, + step=day, + filterType="backwardSmoother", + ) # Forecast the result based on filtered chains def _predictInSample(self, date, days=1): - """ Predict the model's status based on the model of a specific day + """Predict the model's status based on the model of a specific day Args: date: the date the prediction is based on @@ -373,7 +377,7 @@ def _predictInSample(self, date, days=1): """ if date + days > self.n - 1: - raise ValueError('The range is out of sample.') + raise ValueError("The range is out of sample.") predictedObs = [None] * days predictedObsVar = [None] * days @@ -382,8 +386,10 @@ def _predictInSample(self, date, days=1): self.builder.model.prediction.step = 0 for i in range(1, days): # update the evaluation vector - if len(self.builder.dynamicComponents) > 0 or \ - len(self.builder.automaticComponents) > 0: + if ( + len(self.builder.dynamicComponents) > 0 + or len(self.builder.automaticComponents) > 0 + ): self.builder.updateEvaluation(date + i, self.padded_data) self.Filter.predict(self.builder.model) @@ -392,46 +398,38 @@ def _predictInSample(self, date, days=1): return (predictedObs, predictedObsVar) -# =========================== model helper function ========================== - + # =========================== model helper function ========================== # to set model to a specific date def _setModelStatus(self, date=0): - """ Set the model status to a specific date (the date mush have been filtered) - - """ - if date < self.result.filteredSteps[0] or \ - date > self.result.filteredSteps[1]: - raise ValueError('The date has yet to be filtered yet. ' - 'Check the in object.') - - self._reverseCopy(model=self.builder.model, - result=self.result, - step=date) - if len(self.builder.dynamicComponents) > 0 or \ - len(self.builder.automaticComponents) > 0: + """Set the model status to a specific date (the date mush have been filtered)""" + if date < self.result.filteredSteps[0] or date > self.result.filteredSteps[1]: + raise ValueError( + "The date has yet to be filtered yet. " + "Check the in object." + ) + + self._reverseCopy(model=self.builder.model, result=self.result, step=date) + if ( + len(self.builder.dynamicComponents) > 0 + or len(self.builder.automaticComponents) > 0 + ): self.builder.updateEvaluation(date, self.padded_data) - # reset model to initial status def _resetModelStatus(self): - """ Reset the model to the prior status - - """ + """Reset the model to the prior status""" self.builder.model.state = self.builder.statePrior self.builder.model.sysVar = self.builder.sysVarPrior self.builder.model.noiseVar = self.builder.noiseVar self.builder.model.df = self.builder.initialDegreeFreedom self.builder.model.initializeObservation() - # a function used to copy result from the model to the result def _copy(self, model, result, step, filterType): - """ Copy result from the model to _result class - - """ + """Copy result from the model to _result class""" - if filterType == 'forwardFilter': + if filterType == "forwardFilter": result.filteredObs[step] = model.obs result.predictedObs[step] = model.prediction.obs result.filteredObsVar[step] = model.obsVar @@ -446,17 +444,14 @@ def _copy(self, model, result, step, filterType): if self.data[step] is None: self.padded_data[step] = result.filteredObs[step] - elif filterType == 'backwardSmoother': + elif filterType == "backwardSmoother": result.smoothedState[step] = model.state result.smoothedObs[step] = model.obs result.smoothedCov[step] = model.sysVar result.smoothedObsVar[step] = model.obsVar - def _reverseCopy(self, model, result, step): - """ Copy result from _result class to the model - - """ + """Copy result from _result class to the model""" model.obs = result.filteredObs[step] model.prediction.obs = result.predictedObs[step] @@ -469,29 +464,24 @@ def _reverseCopy(self, model, result, step): model.noiseVar = result.noiseVar[step] model.df = result.df[step] - # check if the data size matches the dynamic features def _checkFeatureSize(self): - """ Check features's n matches the data's n - - """ + """Check features's n matches the data's n""" if len(self.builder.dynamicComponents) > 0: for name in self.builder.dynamicComponents: if self.builder.dynamicComponents[name].n != self.n: - raise ValueError(f"The data size of dlm and { name } does not match") - + raise ValueError( + f"The data size of dlm and { name } does not match" + ) def _1DmatrixToArray(self, arrayOf1dMatrix): - """ Change an array of 1 x 1 matrix to normal array. - - """ - to_array = lambda x : None if x is None else x.tolist()[0][0] + """Change an array of 1 x 1 matrix to normal array.""" + to_array = lambda x: None if x is None else x.tolist()[0][0] return [to_array(item) for item in arrayOf1dMatrix] - # function to judge whether a component is in the model def _checkComponent(self, name): - """ Check whether a component is contained by the dlm. + """Check whether a component is contained by the dlm. Args: name: the name of the component @@ -499,17 +489,18 @@ def _checkComponent(self, name): Returns: True or error. """ - if name in self.builder.staticComponents or \ - name in self.builder.dynamicComponents or \ - name in self.builder.automaticComponents: + if ( + name in self.builder.staticComponents + or name in self.builder.dynamicComponents + or name in self.builder.automaticComponents + ): return True else: - raise ValueError('No such component.') - + raise ValueError("No such component.") # function to fetch a component def _fetchComponent(self, name): - """ Get the component if the componeng is in the dlm + """Get the component if the componeng is in the dlm Args: name: the name of the component @@ -524,13 +515,12 @@ def _fetchComponent(self, name): elif name in self.builder.automaticComponents: comp = self.builder.automaticComponents[name] else: - raise ValueError('No such component.') + raise ValueError("No such component.") return comp - # check start and end dates that has been filtered on def _checkAndGetWorkingDates(self, filterType): - """ Check the filter status and return the dates that have + """Check the filter status and return the dates that have been filtered according to the filterType. Args: @@ -540,33 +530,36 @@ def _checkAndGetWorkingDates(self, filterType): Returns: a tuple of (start_date, end_date) """ - if filterType == 'forwardFilter' or filterType == 'predict': + if filterType == "forwardFilter" or filterType == "predict": if self.result.filteredSteps != [0, self.n - 1]: - self._logger.info('The fitlered dates are from ' + - str(self.result.filteredSteps[0]) + - ' to ' + str(self.result.filteredSteps[1])) + self._logger.info( + "The fitlered dates are from " + + str(self.result.filteredSteps[0]) + + " to " + + str(self.result.filteredSteps[1]) + ) start = self.result.filteredSteps[0] end = self.result.filteredSteps[1] - elif filterType == 'backwardSmoother': + elif filterType == "backwardSmoother": if self.result.smoothedSteps != [0, self.n - 1]: - self._logger.info('The smoothed dates are from ' + - str(self.result.smoothedSteps[0]) + - ' to ' + str(self.result.smoothedSteps[1])) + self._logger.info( + "The smoothed dates are from " + + str(self.result.smoothedSteps[0]) + + " to " + + str(self.result.smoothedSteps[1]) + ) start = self.result.smoothedSteps[0] end = self.result.smoothedSteps[1] else: - raise ValueError('Incorrect filter type.') + raise ValueError("Incorrect filter type.") return (start, end) - # check the filter status and automatic turn off some plots def _checkPlotOptions(self): - """ Check the filter status and determine the plot options. - - """ + """Check the filter status and determine the plot options.""" # adjust the option according to the filtering status if self.result.filteredSteps[1] == -1: self.options.plotFilteredData = False @@ -575,23 +568,20 @@ def _checkPlotOptions(self): if self.result.smoothedSteps[1] == -1: self.options.plotSmoothedData = False - - def setLoggingLevel(self, level='INFO'): - """ The level of logs that should be printed. - - """ - if level == 'DEBUG': + def setLoggingLevel(self, level="INFO"): + """The level of logs that should be printed.""" + if level == "DEBUG": self._logger.setLevel(logging.DEBUG) - elif level == 'INFO': + elif level == "INFO": self._logger.setLevel(logging.INFO) - elif level == 'WARNING': + elif level == "WARNING": self._logger.setLevel(logging.WARNING) - elif level == 'CRITICAL': + elif level == "CRITICAL": self._logger.setLevel(logging.CRITICAL) else: raise ValueError(f"Log at {level} is not supported.") def getLoggingLevel(self): - """ Get the level of logger.""" + """Get the level of logger.""" return logging.getLevelName(self._logger.getEffectiveLevel()) diff --git a/pydlm/dlm.py b/pydlm/dlm.py index f8c018ed..6af16bda 100644 --- a/pydlm/dlm.py +++ b/pydlm/dlm.py @@ -40,7 +40,7 @@ class dlm(dlmPlotModule, dlmPredictModule, dlmAccessModule, dlmTuneModule): - """ The main class of the dynamic linear model. + """The main class of the dynamic linear model. This is the main class of the Bayeisan dynamic linear model. @@ -76,7 +76,6 @@ class dlm(dlmPlotModule, dlmPredictModule, dlmAccessModule, dlmTuneModule): """ - # define the basic members # initialize the result def __init__(self, data, **options): @@ -88,23 +87,21 @@ def __init__(self, data, **options): # being changed and behaving abnormally. self._predictModel = None - def exportModel(self): - """ Export the dlm builder. Currently the method only support dlm without + """Export the dlm builder. Currently the method only support dlm without dynamic components. """ if length(self.builder.dynamicComponents) > 0: - raise ValueError('Cannot export dlm builder with dynamic components.') + raise ValueError("Cannot export dlm builder with dynamic components.") if not self.initialized: - raise ValueError('Cannot export dlm before the model was initilized.') + raise ValueError("Cannot export dlm before the model was initilized.") return deepcopy(self.builder) - def buildFromModel(self, model): - """ Construct the dlm with exported model from other DLM with status. + """Construct the dlm with exported model from other DLM with status. Args: model: The exported model from other dlm. Must be the return from @@ -113,12 +110,11 @@ def buildFromModel(self, model): """ self._initializeFromBuilder(exported_builder=model) -# ===================== modeling components ===================== - + # ===================== modeling components ===================== # add component def add(self, component): - """ Add new modeling component to the dlm. + """Add new modeling component to the dlm. Currently support: trend, seasonality, autoregression and dynamic component. @@ -134,24 +130,19 @@ def add(self, component): """ self.__add__(component) - def __add__(self, component): self.builder.__add__(component) self.initialized = False return self - # list all components def ls(self): - """ List out all existing components - - """ + """List out all existing components""" self.builder.ls() - # delete one component def delete(self, name): - """ Delete model component by its name + """Delete model component by its name Args: name: the name of the component. @@ -160,11 +151,10 @@ def delete(self, name): self.builder.delete(name) self.initialized = False -# ========================== model training component ======================= - + # ========================== model training component ======================= def fitForwardFilter(self, useRollingWindow=False, windowLength=3): - """ Fit forward filter on the available data. + """Fit forward filter on the available data. User can choose use rolling windowFront or not. If user choose not to use the rolling window, @@ -185,10 +175,10 @@ def fitForwardFilter(self, useRollingWindow=False, windowLength=3): if not self.initialized: self._initialize() - self._logger.info('Starting forward filtering...') + self._logger.info("Starting forward filtering...") if not useRollingWindow: # we start from the last step of previous filtering - if self.result.filteredType == 'non-rolling': + if self.result.filteredType == "non-rolling": start = self.result.filteredSteps[1] + 1 else: start = 0 @@ -197,12 +187,10 @@ def fitForwardFilter(self, useRollingWindow=False, windowLength=3): self.result.smoothedSteps = [0, -1] # determine whether renew should be used - self._forwardFilter(start=start, - end=self.n - 1, - renew=self.options.stable) - self.result.filteredType = 'non-rolling' + self._forwardFilter(start=start, end=self.n - 1, renew=self.options.stable) + self.result.filteredType = "non-rolling" else: - if self.result.filteredType == 'rolling': + if self.result.filteredType == "rolling": windowFront = self.result.filteredSteps[1] + 1 else: windowFront = 0 @@ -210,29 +198,32 @@ def fitForwardFilter(self, useRollingWindow=False, windowLength=3): # backward smoother as well. self.result.smoothedSteps = [0, -1] - self.result.filteredType = 'rolling' + self.result.filteredType = "rolling" # if end is still within (0, windowLength - 1), we should run the # usual ff from if windowFront < windowLength: - self._forwardFilter(start=self.result.filteredSteps[1] + 1, - end=min(windowLength - 1, self.n - 1)) + self._forwardFilter( + start=self.result.filteredSteps[1] + 1, + end=min(windowLength - 1, self.n - 1), + ) # for the remaining date, we use a rolling window for today in range(max(windowFront, windowLength), self.n): - self._forwardFilter(start=today - windowLength + 1, - end=today, - save=today, - ForgetPrevious=True) + self._forwardFilter( + start=today - windowLength + 1, + end=today, + save=today, + ForgetPrevious=True, + ) self.result.filteredSteps = [0, self.n - 1] - self.turnOn('filtered plot') - self.turnOn('predict plot') - - self._logger.info('Forward filtering completed.') + self.turnOn("filtered plot") + self.turnOn("predict plot") + self._logger.info("Forward filtering completed.") def fitBackwardSmoother(self, backLength=None): - """ Fit backward smoothing on the data. Starting from the last observed date. + """Fit backward smoothing on the data. Starting from the last observed date. Args: backLength: integer, indicating how many days the backward smoother @@ -242,27 +233,33 @@ def fitBackwardSmoother(self, backLength=None): # see if the model has been initialized if not self.initialized: - raise ValueError('Backward Smoother has to be run after' + - ' forward filter') + raise ValueError( + "Backward Smoother has to be run after" + " forward filter" + ) if self.result.filteredSteps[1] != self.n - 1: - raise ValueError('Forward Fiter needs to run on full data before' + - 'using backward Smoother') + raise ValueError( + "Forward Fiter needs to run on full data before" + + "using backward Smoother" + ) # default value for backLength if backLength is None: backLength = self.n - self._logger.info('Starting backward smoothing...') + self._logger.info("Starting backward smoothing...") # if the smoothed dates has already been done, we do nothing - if self.result.smoothedSteps[1] == self.n - 1 and \ - self.result.smoothedSteps[0] <= self.n - 1 - backLength + 1: + if ( + self.result.smoothedSteps[1] == self.n - 1 + and self.result.smoothedSteps[0] <= self.n - 1 - backLength + 1 + ): return None # if the smoothed dates start from n - 1, we just need to continue elif self.result.smoothedSteps[1] == self.n - 1: - self._backwardSmoother(start=self.result.smoothedSteps[0] - 1, - days=backLength) + self._backwardSmoother( + start=self.result.smoothedSteps[0] - 1, days=backLength + ) # if the smoothed dates are even earlier, # we need to start from the beginning @@ -270,25 +267,20 @@ def fitBackwardSmoother(self, backLength=None): self._backwardSmoother(start=self.n - 1, days=backLength) self.result.smoothedSteps = [self.n - backLength, self.n - 1] - self.turnOn('smoothed plot') - - self._logger.info('Backward smoothing completed.') + self.turnOn("smoothed plot") + self._logger.info("Backward smoothing completed.") def fit(self): - """ An easy caller for fitting both the forward filter and backward smoother. - - """ + """An easy caller for fitting both the forward filter and backward smoother.""" self.fitForwardFilter() self.fitBackwardSmoother() - -# ======================= data appending, popping and altering =============== - + # ======================= data appending, popping and altering =============== # Append new data or features to the dlm - def append(self, data, component='main'): - """ Append the new data to the main data or the components (new feature data) + def append(self, data, component="main"): + """Append the new data to the main data or the components (new feature data) Args: data: the new data @@ -301,9 +293,9 @@ def append(self, data, component='main'): # initialize the model to ease the modification if not self.initialized: self._initialize() - + # if we are adding new data to the time series - if component == 'main': + if component == "main": # add the data to the self.data self.data.extend(list(data)) @@ -318,8 +310,10 @@ def append(self, data, component='main'): # give a warning to remind to append dynamic components if len(self.builder.dynamicComponents) > 0: - self._logger.warning('Dynamic components need to manually append new data ' - 'when main data has been appended.') + self._logger.warning( + "Dynamic components need to manually append new data " + "when main data has been appended." + ) # if we are adding new data to the features of dynamic components elif component in self.builder.dynamicComponents: @@ -327,20 +321,18 @@ def append(self, data, component='main'): comp.appendNewData(data) else: - raise ValueError('Such dynamic component does not exist.') - + raise ValueError("Such dynamic component does not exist.") # pop the data of a specific date out def popout(self, date): - """ Pop out the data for a given date + """Pop out the data for a given date Args: date: the index indicates which date to be popped out. """ if date < 0 or date > self.n - 1: - raise NameError('The date should be between 0 and ' + - str(self.n - 1)) + raise NameError("The date should be between 0 and " + str(self.n - 1)) # initialize the model to ease the modification if not self.initialized: @@ -369,10 +361,9 @@ def popout(self, date): elif self.result.smoothedSteps[0] > self.result.smoothedSteps[1]: self.result.smoothedSteps = [0, -1] - # alter the data of a specific days - def alter(self, date, data, component='main'): - """ To alter the data for a specific date and a specific component. + def alter(self, date, data, component="main"): + """To alter the data for a specific date and a specific component. Args: date: the date of the altering data @@ -386,8 +377,7 @@ def alter(self, date, data, component='main'): """ if date < 0 or date > self.n - 1: - raise NameError('The date should be between 0 and ' + - str(self.n - 1)) + raise NameError("The date should be between 0 and " + str(self.n - 1)) # initialize the model to ease the modification if not self.initialized: @@ -396,7 +386,7 @@ def alter(self, date, data, component='main'): # To alter the data for the observed time series. # No need to alter `autoReg` or `longSeason` when only the main data # is altered. - if component == 'main': + if component == "main": self.data[date] = data # to alter the feature of a component @@ -405,7 +395,7 @@ def alter(self, date, data, component='main'): comp.alter(date, data) else: - raise NameError('Such dynamic component does not exist.') + raise NameError("Such dynamic component does not exist.") # update the filtered and the smoothed steps self.result.filteredSteps[1] = date - 1 @@ -418,45 +408,39 @@ def alter(self, date, data, component='main'): elif self.result.smoothedSteps[0] > self.result.smoothedSteps[1]: self.result.smoothedSteps = [0, -1] - # ignore the data of a given date def ignore(self, date): - """ Ignore the data for a specific day. treat it as missing data + """Ignore the data for a specific day. treat it as missing data Args: date: the date to ignore. """ if date < 0 or date > self.n - 1: - raise NameError('The date should be between 0 and ' + - str(self.n - 1)) - - self.alter(date=date, data=None, component='main') + raise NameError("The date should be between 0 and " + str(self.n - 1)) -# ================================ control options ========================= + self.alter(date=date, data=None, component="main") + # ================================ control options ========================= def showOptions(self): - """ Print out all the option values - - """ + """Print out all the option values""" allItems = vars(self.options) for item in allItems: - print(item + ': ' + str(allItems[item])) - + print(item + ": " + str(allItems[item])) def stableMode(self, use=True): - """ Turn on the stable mode, i.e., using the renewal strategy. - - Indicate whether the renew strategy should be used to add numerical - stability. When the filter goes over certain steps, - the information contribution of the previous data has decayed - to minimum. In the stable mode, We then ignore those days and - refit the time series starting from current - renewTerm, where - renewTerm is computed according to the discount. Thus, - the effective sample size of the dlm is twice - renewTerm. When discount = 1, there will be no renewTerm, - since all the information will be passed along. + """Turn on the stable mode, i.e., using the renewal strategy. + + Indicate whether the renew strategy should be used to add numerical + stability. When the filter goes over certain steps, + the information contribution of the previous data has decayed + to minimum. In the stable mode, We then ignore those days and + refit the time series starting from current - renewTerm, where + renewTerm is computed according to the discount. Thus, + the effective sample size of the dlm is twice + renewTerm. When discount = 1, there will be no renewTerm, + since all the information will be passed along. """ # if option changes, reset everything if self.options.stable != use: @@ -467,11 +451,10 @@ def stableMode(self, use=True): elif use is False: self.options.stable = False else: - raise ValueError('Incorrect option input') + raise ValueError("Incorrect option input") - - def evolveMode(self, evoType='dependent'): - """ Control whether different component evolve indpendently. If true, + def evolveMode(self, evoType="dependent"): + """Control whether different component evolve indpendently. If true, then the innovation will only be added on each component but not the correlation between the components, so that for component with discount equals to 1, the smoothed results will always be constant. @@ -486,25 +469,23 @@ def evolveMode(self, evoType='dependent'): a dlm object (for chaining purpose) """ # if option changes, reset everything - if (self.options.innovationType == 'whole' and - evoType == 'independent') or \ - (self.options.innovationType == 'component' and - evoType == 'dependent'): + if (self.options.innovationType == "whole" and evoType == "independent") or ( + self.options.innovationType == "component" and evoType == "dependent" + ): self.initialized = False - if evoType == 'independent': - self.options.innovationType = 'component' - elif evoType == 'dependent': - self.options.innovationType = 'whole' + if evoType == "independent": + self.options.innovationType = "component" + elif evoType == "dependent": + self.options.innovationType = "whole" else: - raise NameError('Incorrect option input') + raise NameError("Incorrect option input") # for chaining return self - def noisePrior(self, prior=0): - """ To set the prior for the observational noise. Calling with empty + """To set the prior for the observational noise. Calling with empty argument will enable the auto noise intializer (currently, the min of 1 and the variance of time series). @@ -515,11 +496,11 @@ def noisePrior(self, prior=0): A dlm object (for chaining purpose) """ if prior > 0: - self.options.noise=prior + self.options.noise = prior self.initialized = False else: self.options.useAutoNoise = True - self.initialized = False + self.initialized = False # for chaining return self diff --git a/pydlm/modeler/autoReg.py b/pydlm/modeler/autoReg.py index b549328d..e9d91d88 100644 --- a/pydlm/modeler/autoReg.py +++ b/pydlm/modeler/autoReg.py @@ -11,6 +11,7 @@ similar to @dynamic. """ + import pydlm.base.tools as tl from .component import component @@ -20,7 +21,7 @@ class autoReg(component): - """ The autoReg class allows user to add an autoregressive component to the dlm. + """The autoReg class allows user to add an autoregressive component to the dlm. This code implements the autoregressive component as a child class of component. Different from the dynamic component, the features in the autoReg is generated from the data, and updated according to the data. @@ -62,15 +63,8 @@ class autoReg(component): """ - - def __init__(self, - degree=2, - discount=0.99, - name='ar2', - w=100, - padding=0): - - self.componentType = 'autoReg' + def __init__(self, degree=2, discount=0.99, name="ar2", w=100, padding=0): + self.componentType = "autoReg" self.d = degree self.name = name self.discount = np.ones(self.d) * discount @@ -87,61 +81,52 @@ def __init__(self, self.createCovPrior(scale=w) self.createMeanPrior() - def createEvaluation(self, step, data): - """ The evaluation matrix for auto regressor. - - """ + """The evaluation matrix for auto regressor.""" if step > len(data): raise NameError("There is no sufficient data for creating autoregressor.") # We pad numbers if the step is too early - self.evaluation = np.array([[self.padding] * (self.d - step) + - list(data[max(0, (step - self.d)) : step])]) - + self.evaluation = np.array( + [ + [self.padding] * (self.d - step) + + list(data[max(0, (step - self.d)) : step]) + ] + ) def createTransition(self): - """ Create the transition matrix. + """Create the transition matrix. For the auto regressor component, the transition matrix is just the identity matrix """ self.transition = np.eye(self.d) - - def createCovPrior(self, cov = None, scale = 1e6): - """ Create the prior covariance matrix for the latent states - - """ + def createCovPrior(self, cov=None, scale=1e6): + """Create the prior covariance matrix for the latent states""" if cov is None: self.covPrior = np.eye(self.d) * scale else: self.covPrior = cov * scale - - def createMeanPrior(self, mean = None, scale = 1): - """ Create the prior latent state - - """ + def createMeanPrior(self, mean=None, scale=1): + """Create the prior latent state""" if mean is None: self.meanPrior = np.zeros((self.d, 1)) * scale else: self.meanPrior = mean * scale - def checkDimensions(self): - """ if user supplies their own covPrior and meanPrior, this can + """if user supplies their own covPrior and meanPrior, this can be used to check if the dimension matches """ tl.checker.checkVectorDimension(self.meanPrior, self.covPrior) - def updateEvaluation(self, step, data): self.createEvaluation(step=step, data=data) - def appendNewData(self, data): - """ AutoReg append new data automatically with the main time series. Nothing + """AutoReg append new data automatically with the main time series. Nothing needs to be done here. """ diff --git a/pydlm/modeler/builder.py b/pydlm/modeler/builder.py index 38e95036..ccc341b4 100644 --- a/pydlm/modeler/builder.py +++ b/pydlm/modeler/builder.py @@ -32,7 +32,7 @@ class builder: - """ The main modeling part of a dynamic linear model. It allows the users to + """The main modeling part of a dynamic linear model. It allows the users to custemize their own model. User can add, delete any components like trend or seasonality to the builder, or view the existing components. Builder will finally assemble all components to a big model for further training @@ -74,7 +74,6 @@ class builder: # create members def __init__(self, logger=None): - # the basic model structure for running kalman filter self.model = None self.initialized = False @@ -125,7 +124,7 @@ def __init__(self, logger=None): # The function that allows the user to add components def add(self, component): - """ Add a new model component to the builder. + """Add a new model component to the builder. Args: component: a model component, any class implements @component class @@ -135,72 +134,70 @@ def add(self, component): self.__add__(component) def __add__(self, component): - if component.componentType == 'dynamic': + if component.componentType == "dynamic": if component.name in self.dynamicComponents: - raise NameError('Please rename the component to a' - + ' different name.') + raise NameError("Please rename the component to a" + " different name.") self.dynamicComponents[component.name] = component - if component.componentType == 'autoReg' \ - or component.componentType == 'longSeason': + if ( + component.componentType == "autoReg" + or component.componentType == "longSeason" + ): if component.name in self.automaticComponents: - raise NameError('Please rename the component to a' - + ' different name.') + raise NameError("Please rename the component to a" + " different name.") self.automaticComponents[component.name] = component - if component.componentType == 'trend' \ - or component.componentType == 'seasonality': + if ( + component.componentType == "trend" + or component.componentType == "seasonality" + ): if component.name in self.staticComponents: - raise NameError('Please rename the component' + - ' to a different name.') + raise NameError("Please rename the component" + " to a different name.") self.staticComponents[component.name] = component # we use seasonality's discount to adjust the renewTerm - if component.componentType == 'seasonality': + if component.componentType == "seasonality": if self.renewDiscount is None: self.renewDiscount = 1.0 - self.renewDiscount = min(self.renewDiscount, - min(component.discount)) + self.renewDiscount = min(self.renewDiscount, min(component.discount)) self.initialized = False return self # print all components to the client def ls(self): - """ List out all the existing components to the model - - """ + """List out all the existing components to the model""" if len(self.staticComponents) > 0: - print('The static components are') + print("The static components are") for name in self.staticComponents: comp = self.staticComponents[name] - print(comp.name + ' (degree = ' + str(comp.d) + ')') - print(' ') + print(comp.name + " (degree = " + str(comp.d) + ")") + print(" ") else: - print('There is no static component.') - print(' ') + print("There is no static component.") + print(" ") if len(self.dynamicComponents) > 0: - print('The dynamic components are') + print("The dynamic components are") for name in self.dynamicComponents: comp = self.dynamicComponents[name] - print(comp.name + ' (dimension = ' + str(comp.d) + ')') - print(' ') + print(comp.name + " (dimension = " + str(comp.d) + ")") + print(" ") else: - print('There is no dynamic component.') - print(' ') + print("There is no dynamic component.") + print(" ") if len(self.automaticComponents) > 0: - print('The automatic components are') + print("The automatic components are") for name in self.automaticComponents: comp = self.automaticComponents[name] - print(comp.name + ' (dimension = ' + str(comp.d) + ')') + print(comp.name + " (dimension = " + str(comp.d) + ")") else: - print('There is no automatic component.') + print("There is no automatic component.") # delete the componet that pointed out by the client def delete(self, name): - """ Delete a specific component from dlm by its name. + """Delete a specific component from dlm by its name. Args: name: the name of the component. Can be read from ls() @@ -214,7 +211,7 @@ def delete(self, name): elif name in self.automaticComponents: del self.automaticComponents[name] else: - raise NameError('Such component does not exisit!') + raise NameError("Such component does not exisit!") self.initialized = False @@ -222,21 +219,21 @@ def delete(self, name): # noise is the prior guess of the variance of the observed data # data is used by auto regressor. def initialize(self, data=[], noise=1): - """ Initialize the model. It construct the baseModel by assembling all + """Initialize the model. It construct the baseModel by assembling all quantities from the components. Args: noise: the initial guess of the variance of the observation noise. """ - if len(self.staticComponents) == 0 and \ - len(self.dynamicComponents) == 0 and \ - len(self.automaticComponents) == 0: - - raise NameError('The model must contain at least' + - ' one component') + if ( + len(self.staticComponents) == 0 + and len(self.dynamicComponents) == 0 + and len(self.automaticComponents) == 0 + ): + raise NameError("The model must contain at least" + " one component") # construct transition, evaluation, prior state, prior covariance - self._logger.info('Initializing models...') + self._logger.info("Initializing models...") transition = None evaluation = None state = None @@ -250,8 +247,7 @@ def initialize(self, data=[], noise=1): for i in self.staticComponents: comp = self.staticComponents[i] transition = mt.matrixAddInDiag(transition, comp.transition) - evaluation = mt.matrixAddByCol(evaluation, - comp.evaluation) + evaluation = mt.matrixAddByCol(evaluation, comp.evaluation) state = mt.matrixAddByRow(state, comp.meanPrior) sysVar = mt.matrixAddInDiag(sysVar, comp.covPrior) self.discount = np.concatenate((self.discount, comp.discount)) @@ -265,13 +261,11 @@ def initialize(self, data=[], noise=1): comp = self.dynamicComponents[i] comp.updateEvaluation(0) transition = mt.matrixAddInDiag(transition, comp.transition) - evaluation = mt.matrixAddByCol(evaluation, - comp.evaluation) + evaluation = mt.matrixAddByCol(evaluation, comp.evaluation) state = mt.matrixAddByRow(state, comp.meanPrior) sysVar = mt.matrixAddInDiag(sysVar, comp.covPrior) self.discount = np.concatenate((self.discount, comp.discount)) - self.componentIndex[i] = (currentIndex, - currentIndex + comp.d - 1) + self.componentIndex[i] = (currentIndex, currentIndex + comp.d - 1) currentIndex += comp.d # if the model contains the automatic dynamic part, we add @@ -282,24 +276,24 @@ def initialize(self, data=[], noise=1): comp = self.automaticComponents[i] comp.updateEvaluation(0, data) transition = mt.matrixAddInDiag(transition, comp.transition) - evaluation = mt.matrixAddByCol(evaluation, - comp.evaluation) + evaluation = mt.matrixAddByCol(evaluation, comp.evaluation) state = mt.matrixAddByRow(state, comp.meanPrior) sysVar = mt.matrixAddInDiag(sysVar, comp.covPrior) self.discount = np.concatenate((self.discount, comp.discount)) - self.componentIndex[i] = (currentIndex, - currentIndex + comp.d - 1) + self.componentIndex[i] = (currentIndex, currentIndex + comp.d - 1) currentIndex += comp.d self.statePrior = state self.sysVarPrior = sysVar self.noiseVar = np.array([[noise]]) - self.model = baseModel(transition=transition, - evaluation=evaluation, - noiseVar=np.array([[noise]]), - sysVar=sysVar, - state=state, - df=self.initialDegreeFreedom) + self.model = baseModel( + transition=transition, + evaluation=evaluation, + noiseVar=np.array([[noise]]), + sysVar=sysVar, + state=state, + df=self.initialDegreeFreedom, + ) self.model.initializeObservation() # compute the renew period @@ -307,11 +301,12 @@ def initialize(self, data=[], noise=1): self.renewDiscount = np.min(self.discount) if self.renewDiscount < 1.0 - 1e-8: - self.renewTerm = np.log(0.001 * (1 - self.renewDiscount)) \ - / np.log(self.renewDiscount) + self.renewTerm = np.log(0.001 * (1 - self.renewDiscount)) / np.log( + self.renewDiscount + ) self.initialized = True - self._logger.info('Initialization finished.') + self._logger.info("Initialization finished.") # Initialize from another builder exported from other dlm class def initializeFromBuilder(self, data, exported_builder): @@ -338,18 +333,19 @@ def initializeFromBuilder(self, data, exported_builder): self.renewDiscount = np.min(self.discount) if self.renewDiscount < 1.0 - 1e-8: - self.renewTerm = np.log(0.001 * (1 - self.renewDiscount)) \ - / np.log(self.renewDiscount) - + self.renewTerm = np.log(0.001 * (1 - self.renewDiscount)) / np.log( + self.renewDiscount + ) + self.initialized = True - self._logger.info('Initialization finished.') - + self._logger.info("Initialization finished.") + # This function allows the model to update the dynamic evaluation vector, # so that the model can handle control variables # This function should be called only when dynamicComponents is not empty # data is used by auto regressor. def updateEvaluation(self, step, data): - """ Update the evaluation matrix of the model to a specific date. + """Update the evaluation matrix of the model to a specific date. It loops over all dynamic components and update their evaluation matrix and then reconstruct the model evaluation matrix by incorporating the new evaluations @@ -369,11 +365,13 @@ def updateEvaluation(self, step, data): for i in self.dynamicComponents: comp = self.dynamicComponents[i] comp.updateEvaluation(step) - self.model.evaluation[0, self.componentIndex[i][0]: - (self.componentIndex[i][1] + 1)] = comp.evaluation + self.model.evaluation[ + 0, self.componentIndex[i][0] : (self.componentIndex[i][1] + 1) + ] = comp.evaluation for i in self.automaticComponents: comp = self.automaticComponents[i] comp.updateEvaluation(step, data) - self.model.evaluation[0, self.componentIndex[i][0]: - (self.componentIndex[i][1] + 1)] = comp.evaluation + self.model.evaluation[ + 0, self.componentIndex[i][0] : (self.componentIndex[i][1] + 1) + ] = comp.evaluation diff --git a/pydlm/modeler/component.py b/pydlm/modeler/component.py index c4973a94..b0cec1b3 100644 --- a/pydlm/modeler/component.py +++ b/pydlm/modeler/component.py @@ -10,25 +10,29 @@ components based on this abstract class. """ + from abc import ABCMeta, abstractmethod +import numpy as np + # We define an abstract class which can further be used # to create different types of model components, inclusing # trend, seasonality and other structures class component: - """ The abstract class provides the basic structure for all model components - + """The abstract class provides the basic structure for all model components + Methods: createEvaluation: create the initial evaluation matrix createTransition: create the initial transition matrix createCovPrior: create a simple prior covariance matrix createMeanPrior: create a simple prior latent state - checkDimensions: if user supplies their own covPrior and meanPrior, this can + checkDimensions: if user supplies their own covPrior and meanPrior, this can be used to check if the dimension matches - + """ + __metaclass__ = ABCMeta def __init__(self): @@ -47,59 +51,58 @@ def __eq__(self, other): if not isinstance(other, component): return NotImplemented else: - return (self.equalOrNone(self.d, other.d) and - self.equalOrNone(self.name, other.name) and - self.equalOrNone(self.componentType, other.componentType) and - self.npEqualOrNone(self.discount, other.discount) and - self.npEqualOrNone(self.evaluation, other.evaluation) and - self.npEqualOrNone(self.transition, other.transition) and - self.npEqualOrNone(self.covPrior, other.covPrior) and - self.npEqualOrNone(self.meanPrior, other.meanPrior)) - + return ( + np.array_equal(self.d, other.d) + and np.array_equal(self.name, other.name) + and np.array_equal(self.componentType, other.componentType) + and np.array_equal(self.discount, other.discount) + and np.array_equal(self.evaluation, other.evaluation) + and np.array_equal(self.transition, other.transition) + and np.array_equal(self.covPrior, other.covPrior) + and np.array_equal(self.meanPrior, other.meanPrior) + ) # define the evaluation matrix for the component @abstractmethod - def createEvaluation(self): pass + def createEvaluation(self): + pass + """ Create the evaluation matrix """ - # define the transition matrix for the component @abstractmethod - def createTransition(self): pass + def createTransition(self): + pass + """ Create the transition matrix """ - # define the prior distribution for the covariance for the component @abstractmethod - def createCovPrior(self): pass + def createCovPrior(self): + pass + """ Create the prior covariance matrix for the latent states """ - # define the prior distribution for the mean vector for the component @abstractmethod - def createMeanPrior(self): pass + def createMeanPrior(self): + pass + """ Create the prior latent state """ - # check the matrix dimensions in case user supplied matrices are wrong @abstractmethod - def checkDimensions(self): pass + def checkDimensions(self): + pass + """ Check the dimensionality of the state and covariance """ - - def equalOrNone(self, a, b): - """Check if a and b are equal or both are None""" - return (a is None and b is None) or a == b - - def npEqualOrNone(self, a, b): - """Check if a and b are equal or both are None for NP arrays""" - return (a is None and b is None) or (a == b).all() diff --git a/pydlm/modeler/dynamic.py b/pydlm/modeler/dynamic.py index 16ec98cd..435e1f6f 100644 --- a/pydlm/modeler/dynamic.py +++ b/pydlm/modeler/dynamic.py @@ -14,6 +14,7 @@ The name dynamic means that the features are changing over time. """ + import numpy as np from collections.abc import MutableSequence from copy import deepcopy @@ -25,7 +26,7 @@ class dynamic(component): - """ The dynamic component that allows user to add controlled variables, + """The dynamic component that allows user to add controlled variables, providing one building block for the dynamic linear model. It decribes a dynamic component in the time series data. It basically allows user to supply covariate or controlled variable into the dlm, @@ -64,24 +65,19 @@ class dynamic(component): """ - - def __init__(self, - features = None, - discount = 0.99, - name = 'dynamic', - w=100): - + def __init__(self, features=None, discount=0.99, name="dynamic", w=100): self.n = len(features) self.d = len(features[0]) if self.hasMissingData(features): - raise ValueError("The current version does not support missing data" + - "in the features.") + raise ValueError( + "The current version does not support missing data" + "in the features." + ) self.features = deepcopy(features) if isinstance(features, np.matrix) or isinstance(features, np.ndarray): self.features = self.features.tolist() - self.componentType = 'dynamic' + self.componentType = "dynamic" self.name = name self.discount = np.ones(self.d) * discount @@ -97,47 +93,38 @@ def __init__(self, self.createCovPrior(scale=w) self.createMeanPrior() - def createEvaluation(self, step): - """ The evaluation matrix for the dynamic component change over time. + """The evaluation matrix for the dynamic component change over time. It equals to the value of the features or the controlled variables at a given date """ self.evaluation = np.array([self.features[step]]) - def createTransition(self): - """ Create the transition matrix. + """Create the transition matrix. For the dynamic component, the transition matrix is just the identity matrix """ self.transition = np.eye(self.d) - - def createCovPrior(self, cov = None, scale = 1e6): - """ Create the prior covariance matrix for the latent states - - """ + def createCovPrior(self, cov=None, scale=1e6): + """Create the prior covariance matrix for the latent states""" if cov is None: self.covPrior = np.eye(self.d) * scale else: self.covPrior = cov * scale - - def createMeanPrior(self, mean = None, scale = 1): - """ Create the prior latent state - - """ + def createMeanPrior(self, mean=None, scale=1): + """Create the prior latent state""" if mean is None: self.meanPrior = np.zeros((self.d, 1)) * scale else: self.meanPrior = mean * scale - def checkDimensions(self): - """ if user supplies their own covPrior and meanPrior, this can + """if user supplies their own covPrior and meanPrior, this can be used to check if the dimension matches """ @@ -146,9 +133,7 @@ def checkDimensions(self): # Recursively heck if there is any none data. We currently don't support # missing data for features. def hasMissingData(self, features): - """ Check whether the list contains None - - """ + """Check whether the list contains None""" for item in features: if isinstance(item, MutableSequence): if self.hasMissingData(item): @@ -158,9 +143,8 @@ def hasMissingData(self, features): return True return False - def updateEvaluation(self, step): - """ update the evaluation matrix to a specific date + """update the evaluation matrix to a specific date This function is used when fitting the forward filter and backward smoother in need of updating the correct evaluation matrix @@ -169,11 +153,10 @@ def updateEvaluation(self, step): self.evaluation = np.array([self.features[step]]) self.step = step else: - raise ValueError('The step is out of range') - + raise ValueError("The step is out of range") def appendNewData(self, newData): - """ For updating feature matrix when new data is added. + """For updating feature matrix when new data is added. Args: newData: is a list of list. The inner list is the feature vector. The outer @@ -181,15 +164,15 @@ def appendNewData(self, newData): """ if self.hasMissingData(newData): - raise ValueError("The current version does not support missing data" + - "in the features.") - + raise ValueError( + "The current version does not support missing data" + "in the features." + ) + self.features.extend(tl.duplicateList(newData)) self.n = len(self.features) - def popout(self, date): - """ For deleting the feature data of a specific date. + """For deleting the feature data of a specific date. Args: date: the index of which to be deleted. @@ -198,9 +181,8 @@ def popout(self, date): self.features.pop(date) self.n -= 1 - def alter(self, date, feature): - """ Change the corresponding + """Change the corresponding feature matrix. Args: @@ -209,8 +191,8 @@ def alter(self, date, feature): """ if self.hasMissingData(feature): - raise ValueError("The current version does not support missing data" + - "in the features.") + raise ValueError( + "The current version does not support missing data" + "in the features." + ) else: self.features[date] = feature - diff --git a/pydlm/modeler/longSeason.py b/pydlm/modeler/longSeason.py index 501e976f..db854507 100644 --- a/pydlm/modeler/longSeason.py +++ b/pydlm/modeler/longSeason.py @@ -21,13 +21,14 @@ similar to @dynamic. """ + from .autoReg import autoReg import logging import numpy as np class longSeason(autoReg): - """ The longSeason class alows user to add a long seasonality component + """The longSeason class alows user to add a long seasonality component to the dlm. The difference between the long seasonality is that 1) The seasonality component use each date as a unit and change in a given periodicity. For example, 1, 2, 3, 4, 1, 2, 3, 4. @@ -61,31 +62,20 @@ class longSeason(autoReg): """ - - def __init__(self, - period=4, - stay=7, - discount=0.99, - name='longSeason', - w=100): - + def __init__(self, period=4, stay=7, discount=0.99, name="longSeason", w=100): self.period = period self.stay = stay - super().__init__(degree=period, - discount=discount, - name=name, - w=w) + super().__init__(degree=period, discount=discount, name=name, w=w) # modify the type to be longSeason - self.componentType = 'longSeason' + self.componentType = "longSeason" # Initialize the evaluation vector self.evaluation = np.array([[0] * self.period]) - def updateEvaluation(self, step, data=None): - """ update the evaluation matrix to a specific date + """update the evaluation matrix to a specific date This function is used when fitting the forward filter and backward smoother in need of updating the correct evaluation matrix diff --git a/pydlm/modeler/matrixTools.py b/pydlm/modeler/matrixTools.py index 03de46d4..8960aec3 100644 --- a/pydlm/modeler/matrixTools.py +++ b/pydlm/modeler/matrixTools.py @@ -17,7 +17,6 @@ def matrixAddInDiag(A, B): newMatrixB = np.concatenate((np.zeros((Bn, Ap)), B), 1) return np.concatenate((newMatrixA, newMatrixB), 0) - # A + B = (A; B) @staticmethod def matrixAddByRow(A, B): @@ -28,7 +27,6 @@ def matrixAddByRow(A, B): else: return np.concatenate((A, B), 0) - # A + B = (A B) @staticmethod def matrixAddByCol(A, B): @@ -39,7 +37,6 @@ def matrixAddByCol(A, B): else: return np.concatenate((A, B), 1) - @staticmethod def AddTwoVectors(a, b): if a is None: diff --git a/pydlm/modeler/seasonality.py b/pydlm/modeler/seasonality.py index a64402b9..46322ac3 100644 --- a/pydlm/modeler/seasonality.py +++ b/pydlm/modeler/seasonality.py @@ -13,6 +13,7 @@ or cos relationship between each state. They can be arbitrarily valued. """ + import numpy as np from .component import component import pydlm.base.tools as tl @@ -60,17 +61,11 @@ class seasonality(component): """ - - def __init__(self, - period = 7, - discount = 0.99, - name = 'seasonality', - w=100): - + def __init__(self, period=7, discount=0.99, name="seasonality", w=100): if period <= 1: - raise ValueError('Period has to be greater than 1.') + raise ValueError("Period has to be greater than 1.") self.d = period - self.componentType = 'seasonality' + self.componentType = "seasonality" self.name = name self.discount = np.ones(self.d) * discount @@ -89,15 +84,11 @@ def __init__(self, # create form free seasonality component self.freeForm() - def createEvaluation(self): - """ Create the evaluation matrix - - """ + """Create the evaluation matrix""" self.evaluation = np.zeros((1, self.d)) self.evaluation[0, 0] = 1 - # The transition matrix takes special form as # G = [0 1 0] # [0 0 1] @@ -105,7 +96,7 @@ def createEvaluation(self): # So everyt time, when G * G, we rotate the vector once, which results # in the seasonality performance def createTransition(self): - """ Create the transition matrix. + """Create the transition matrix. According to Hurrison and West (1999), the transition matrix of seasonality takes a form of\n @@ -119,42 +110,34 @@ def createTransition(self): self.transition = np.diag(np.ones(self.d - 1), 1) self.transition[self.d - 1, 0] = 1 - - def createCovPrior(self, cov = 1e7): - """Create the prior covariance matrix for the latent states. - - """ + def createCovPrior(self, cov=1e7): + """Create the prior covariance matrix for the latent states.""" self.covPrior = np.eye(self.d) * cov - - def createMeanPrior(self, mean = 0): - """ Create the prior latent state - - """ + def createMeanPrior(self, mean=0): + """Create the prior latent state""" self.meanPrior = np.ones((self.d, 1)) * mean - # Form free seasonality component ensures that sum(mean) = 0 # We use the formular from "Bayesian forecasting and dynamic linear models" # Page 242 def freeForm(self): - """ The technique used in Hurrison and West (1999). After calling this method, + """The technique used in Hurrison and West (1999). After calling this method, The latent states sum up to 0 and the covariance matrix is degenerate to have rank d - 1, so that the sum of the latent states will never change when the system evolves """ if self.covPrior is None or self.meanPrior is None: - raise NameError('freeForm can only be called after prior created.') + raise NameError("freeForm can only be called after prior created.") else: u = np.sum(self.covPrior) A = np.sum(self.covPrior, axis=1, keepdims=True) / u self.meanPrior = self.meanPrior - A * np.sum(self.meanPrior) self.covPrior = self.covPrior - np.dot(A, A.T) * u - def checkDimensions(self): - """ if user supplies their own covPrior and meanPrior, this can + """if user supplies their own covPrior and meanPrior, this can be used to check if the dimension matches """ diff --git a/pydlm/modeler/trends.py b/pydlm/modeler/trends.py index 6284d723..6323890e 100644 --- a/pydlm/modeler/trends.py +++ b/pydlm/modeler/trends.py @@ -9,6 +9,7 @@ It decribes a latent polynomial trending in the time series data. """ + import numpy as np from .component import component import pydlm.base.tools as tl @@ -18,7 +19,7 @@ class trend(component): - """ The trend component that features the polynomial trending, + """The trend component that features the polynomial trending, providing one building block for the dynamic linear model. It decribes a latent polynomial trending in the time series data. @@ -50,18 +51,12 @@ class trend(component): """ - - def __init__(self, - degree = 0, - discount = 0.99, - name = 'trend', - w=100): - + def __init__(self, degree=0, discount=0.99, name="trend", w=100): if degree < 0: - raise NameError('degree has to be non-negative') + raise NameError("degree has to be non-negative") self.d = degree + 1 self.name = name - self.componentType = 'trend' + self.componentType = "trend" self.discount = np.ones(self.d) * discount # Initialize all basic quantities @@ -76,15 +71,11 @@ def __init__(self, self.createCovPrior(cov=w) self.createMeanPrior() - def createEvaluation(self): - """ Create the evaluation matrix - - """ + """Create the evaluation matrix""" self.evaluation = np.zeros((1, self.d)) self.evaluation[0, 0] = 1 - def createTransition(self): """Create the transition matrix @@ -100,23 +91,16 @@ def createTransition(self): self.transition = np.zeros((self.d, self.d)) self.transition[np.triu_indices(self.d)] = 1 - def createCovPrior(self, cov=1e7): - """Create the prior covariance matrix for the latent states. - - """ + """Create the prior covariance matrix for the latent states.""" self.covPrior = np.eye(self.d) * cov - def createMeanPrior(self, mean=0): - """ Create the prior latent state - - """ + """Create the prior latent state""" self.meanPrior = np.ones((self.d, 1)) * mean - def checkDimensions(self): - """ if user supplies their own covPrior and meanPrior, this can + """if user supplies their own covPrior and meanPrior, this can be used to check if the dimension matches """ diff --git a/pydlm/plot/__init__.py b/pydlm/plot/__init__.py index f741fa7d..8783eed4 100644 --- a/pydlm/plot/__init__.py +++ b/pydlm/plot/__init__.py @@ -1,3 +1,3 @@ -#__all__ = ['dlmPlot'] +# __all__ = ['dlmPlot'] -#import pydlm.plot.dlmPlot as dlmPlot +# import pydlm.plot.dlmPlot as dlmPlot diff --git a/pydlm/plot/dlmPlot.py b/pydlm/plot/dlmPlot.py index 020172d6..0e43702c 100644 --- a/pydlm/plot/dlmPlot.py +++ b/pydlm/plot/dlmPlot.py @@ -8,6 +8,7 @@ This module provide the plotting functionality for dlm """ + import matplotlib.pyplot as plt from pydlm.base.tools import getInterval @@ -25,71 +26,96 @@ def plotInOneFigure(time, data, result, options): options: options for the plot, for details please refer to @dlm """ # plot the original data - plotData(time=time, data=data, - showDataPoint=options.showDataPoint, color=options.dataColor, - label='time series') + plotData( + time=time, + data=data, + showDataPoint=options.showDataPoint, + color=options.dataColor, + label="time series", + ) # plot fitered results if needed if options.plotFilteredData: start = result.filteredSteps[0] end = result.filteredSteps[1] + 1 - plotData(time=time[start:end], - data=to1dArray(result.filteredObs[start:end]), - showDataPoint=options.showFittedPoint, - color=options.filteredColor, - label='filtered series') + plotData( + time=time[start:end], + data=to1dArray(result.filteredObs[start:end]), + showDataPoint=options.showFittedPoint, + color=options.filteredColor, + label="filtered series", + ) if options.showConfidenceInterval: - upper, lower = getInterval(result.filteredObs[start:end], - result.filteredObsVar[start:end], - p=options.confidence) - - plotInterval(time=time[start:end], - upper=to1dArray(upper), lower=to1dArray(lower), - intervalType=options.intervalType, - color=options.filteredColor) + upper, lower = getInterval( + result.filteredObs[start:end], + result.filteredObsVar[start:end], + p=options.confidence, + ) + + plotInterval( + time=time[start:end], + upper=to1dArray(upper), + lower=to1dArray(lower), + intervalType=options.intervalType, + color=options.filteredColor, + ) # plot predicted results if needed if options.plotPredictedData: start = result.filteredSteps[0] end = result.filteredSteps[1] + 1 - plotData(time=time[start:end], - data=to1dArray(result.predictedObs), - showDataPoint=options.showFittedPoint, - color=options.predictedColor, - label='one-day prediction') + plotData( + time=time[start:end], + data=to1dArray(result.predictedObs), + showDataPoint=options.showFittedPoint, + color=options.predictedColor, + label="one-day prediction", + ) if options.showConfidenceInterval: - upper, lower = getInterval(result.predictedObs[start:end], - result.filteredObsVar[start:end], - p=options.confidence) - - plotInterval(time=time[start:end], - upper=to1dArray(upper), lower=to1dArray(lower), - intervalType=options.intervalType, - color=options.predictedColor) + upper, lower = getInterval( + result.predictedObs[start:end], + result.filteredObsVar[start:end], + p=options.confidence, + ) + + plotInterval( + time=time[start:end], + upper=to1dArray(upper), + lower=to1dArray(lower), + intervalType=options.intervalType, + color=options.predictedColor, + ) # plot smoothed results if needed if options.plotSmoothedData: start = result.smoothedSteps[0] end = result.smoothedSteps[1] + 1 - plotData(time=time[start:end], - data=to1dArray(result.smoothedObs), - showDataPoint=options.showFittedPoint, - color=options.smoothedColor, - label='smoothed series') + plotData( + time=time[start:end], + data=to1dArray(result.smoothedObs), + showDataPoint=options.showFittedPoint, + color=options.smoothedColor, + label="smoothed series", + ) if options.showConfidenceInterval: - upper, lower = getInterval(result.smoothedObs[start:end], - result.smoothedObsVar[start:end], - p=options.confidence) + upper, lower = getInterval( + result.smoothedObs[start:end], + result.smoothedObsVar[start:end], + p=options.confidence, + ) - plotInterval(time=time[start:end], - upper=to1dArray(upper), lower=to1dArray(lower), - intervalType=options.intervalType, - color=options.smoothedColor) + plotInterval( + time=time[start:end], + upper=to1dArray(upper), + lower=to1dArray(lower), + intervalType=options.intervalType, + color=options.smoothedColor, + ) - plt.legend(loc='best', shadow=True) # , fontsize = 'x-large') + plt.legend(loc="best", shadow=True) # , fontsize = 'x-large') def plotInMultipleFigure(time, data, result, options): @@ -104,8 +130,9 @@ def plotInMultipleFigure(time, data, result, options): options: options for the plot, for details please refer to @dlm """ # first compute how many plots are needed - numOfPlots = options.plotFilteredData + options.plotPredictedData + \ - options.plotSmoothedData + numOfPlots = ( + options.plotFilteredData + options.plotPredictedData + options.plotSmoothedData + ) size = (numOfPlots, 1) location = 1 @@ -116,29 +143,40 @@ def plotInMultipleFigure(time, data, result, options): subplot(size, location) # plot original data - plotData(time=time, data=data, - showDataPoint=options.showDataPoint, color=options.dataColor, - label='time series') + plotData( + time=time, + data=data, + showDataPoint=options.showDataPoint, + color=options.dataColor, + label="time series", + ) start = result.filteredSteps[0] end = result.filteredSteps[1] + 1 if start < end: - plotData(time=time[start:end], - data=to1dArray(result.filteredObs[start:end]), - showDataPoint=options.showFittedPoint, - color=options.filteredColor, - label='filtered series') + plotData( + time=time[start:end], + data=to1dArray(result.filteredObs[start:end]), + showDataPoint=options.showFittedPoint, + color=options.filteredColor, + label="filtered series", + ) if options.showConfidenceInterval: - upper, lower = getInterval(result.filteredObs[start:end], - result.filteredObsVar[start:end], - p=options.confidence) - - plotInterval(time=time[start:end], - upper=to1dArray(upper), lower=to1dArray(lower), - intervalType=options.intervalType, - color=options.filteredColor) - plt.legend(loc='best', shadow=True) # , fontsize = 'x-large') + upper, lower = getInterval( + result.filteredObs[start:end], + result.filteredObsVar[start:end], + p=options.confidence, + ) + + plotInterval( + time=time[start:end], + upper=to1dArray(upper), + lower=to1dArray(lower), + intervalType=options.intervalType, + color=options.filteredColor, + ) + plt.legend(loc="best", shadow=True) # , fontsize = 'x-large') location += 1 # plot predicted results if needed @@ -147,28 +185,39 @@ def plotInMultipleFigure(time, data, result, options): subplot(size, location) # plot original data - plotData(time=time, data=data, - showDataPoint=options.showDataPoint, - color=options.dataColor, label='time series') + plotData( + time=time, + data=data, + showDataPoint=options.showDataPoint, + color=options.dataColor, + label="time series", + ) start = result.filteredSteps[0] end = result.filteredSteps[1] + 1 if start < end: - plotData(time=time[start:end], - data=to1dArray(result.predictedObs), - showDataPoint=options.showFittedPoint, - color=options.predictedColor, - label='one-day prediction') + plotData( + time=time[start:end], + data=to1dArray(result.predictedObs), + showDataPoint=options.showFittedPoint, + color=options.predictedColor, + label="one-day prediction", + ) if options.showConfidenceInterval: - upper, lower = getInterval(result.predictedObs[start:end], - result.filteredObsVar[start:end], - p=options.confidence) - plotInterval(time=time[start:end], - upper=to1dArray(upper), lower=to1dArray(lower), - intervalType=options.intervalType, - color=options.predictedColor) - plt.legend(loc='best', shadow=True) + upper, lower = getInterval( + result.predictedObs[start:end], + result.filteredObsVar[start:end], + p=options.confidence, + ) + plotInterval( + time=time[start:end], + upper=to1dArray(upper), + lower=to1dArray(lower), + intervalType=options.intervalType, + color=options.predictedColor, + ) + plt.legend(loc="best", shadow=True) location += 1 # plot smoothed results if needed @@ -177,28 +226,40 @@ def plotInMultipleFigure(time, data, result, options): subplot(size, location) # plot original data - plotData(time=time, data=data, - showDataPoint=options.showDataPoint, - color=options.dataColor, label='time series') + plotData( + time=time, + data=data, + showDataPoint=options.showDataPoint, + color=options.dataColor, + label="time series", + ) start = result.smoothedSteps[0] end = result.smoothedSteps[1] + 1 if start < end: - plotData(time=time[start:end], - data=to1dArray(result.smoothedObs), - showDataPoint=options.showFittedPoint, - color=options.smoothedColor, - label='smoothed series') + plotData( + time=time[start:end], + data=to1dArray(result.smoothedObs), + showDataPoint=options.showFittedPoint, + color=options.smoothedColor, + label="smoothed series", + ) if options.showConfidenceInterval: - upper, lower = getInterval(result.smoothedObs[start:end], - result.smoothedObsVar[start:end], - p=options.confidence) - plotInterval(time=time[start:end], - upper=to1dArray(upper), lower=to1dArray(lower), - intervalType=options.intervalType, - color=options.smoothedColor) - plt.legend(loc='best', shadow=True) + upper, lower = getInterval( + result.smoothedObs[start:end], + result.smoothedObsVar[start:end], + p=options.confidence, + ) + plotInterval( + time=time[start:end], + upper=to1dArray(upper), + lower=to1dArray(lower), + intervalType=options.intervalType, + color=options.smoothedColor, + ) + plt.legend(loc="best", shadow=True) + # ============================ plot for latents ============================ @@ -222,8 +283,7 @@ def plotLatentState(time, coordinates, result, options, name): for i, dim in enumerate(coordinates): subplot(size, i + 1) plotSingleState(time, dim, result, options) - plt.title('Filter result for dimension ' + str(i) + - ' in component: ' + name) + plt.title("Filter result for dimension " + str(i) + " in component: " + name) def plotSingleState(time, dimension, result, options): @@ -242,69 +302,87 @@ def plotSingleState(time, dimension, result, options): start = result.filteredSteps[0] end = result.filteredSteps[1] + 1 data = [item[dimension, 0] for item in result.filteredState[start:end]] - var = [abs(item[dimension, dimension]) - for item in result.filteredCov[start:end]] - - plotData(time=time[start:end], - data=data, - showDataPoint=options.showFittedPoint, - color=options.filteredColor, - label='filtered state') + var = [ + abs(item[dimension, dimension]) for item in result.filteredCov[start:end] + ] + + plotData( + time=time[start:end], + data=data, + showDataPoint=options.showFittedPoint, + color=options.filteredColor, + label="filtered state", + ) if options.showConfidenceInterval: upper, lower = getInterval(data, var, p=options.confidence) - plotInterval(time=time[start:end], - upper=upper, lower=lower, - intervalType=options.intervalType, - color=options.filteredColor) + plotInterval( + time=time[start:end], + upper=upper, + lower=lower, + intervalType=options.intervalType, + color=options.filteredColor, + ) # plot predicted results if needed if options.plotPredictedData: start = result.filteredSteps[0] end = result.filteredSteps[1] + 1 - data = [item[dimension, 0] - for item in result.predictedState[start:end]] - var = [abs(item[dimension, dimension]) - for item in result.predictedCov[start:end]] - - plotData(time=time[start:end], - data=data, - showDataPoint=options.showFittedPoint, - color=options.predictedColor, - label='predicted state') + data = [item[dimension, 0] for item in result.predictedState[start:end]] + var = [ + abs(item[dimension, dimension]) for item in result.predictedCov[start:end] + ] + + plotData( + time=time[start:end], + data=data, + showDataPoint=options.showFittedPoint, + color=options.predictedColor, + label="predicted state", + ) if options.showConfidenceInterval: upper, lower = getInterval(data, var, p=options.confidence) - plotInterval(time=time[start:end], - upper=upper, lower=lower, - intervalType=options.intervalType, - color=options.predictedColor) + plotInterval( + time=time[start:end], + upper=upper, + lower=lower, + intervalType=options.intervalType, + color=options.predictedColor, + ) # plot smoothed results if needed if options.plotSmoothedData: start = result.smoothedSteps[0] end = result.smoothedSteps[1] + 1 data = [item[dimension, 0] for item in result.smoothedState[start:end]] - var = [abs(item[dimension, dimension]) - for item in result.smoothedCov[start:end]] - - plotData(time=time[start:end], - data=data, - showDataPoint=options.showFittedPoint, - color=options.smoothedColor, - label='smoothed state') + var = [ + abs(item[dimension, dimension]) for item in result.smoothedCov[start:end] + ] + + plotData( + time=time[start:end], + data=data, + showDataPoint=options.showFittedPoint, + color=options.smoothedColor, + label="smoothed state", + ) if options.showConfidenceInterval: upper, lower = getInterval(data, var, p=options.confidence) - plotInterval(time=time[start:end], - upper=upper, lower=lower, - intervalType=options.intervalType, - color=options.smoothedColor) + plotInterval( + time=time[start:end], + upper=upper, + lower=lower, + intervalType=options.intervalType, + color=options.smoothedColor, + ) + + plt.legend(loc="best", shadow=True) - plt.legend(loc='best', shadow=True) # ============================ plot for component ============================= @@ -329,8 +407,11 @@ def plotComponent(time, data, result, options): """ if options.separatePlot: - numOfPlots = options.plotFilteredData + options.plotPredictedData + \ - options.plotSmoothedData + numOfPlots = ( + options.plotFilteredData + + options.plotPredictedData + + options.plotSmoothedData + ) size = (numOfPlots, 1) location = 1 @@ -342,23 +423,30 @@ def plotComponent(time, data, result, options): start = result.filteredSteps[0] end = result.filteredSteps[1] + 1 - plotData(time=time[start:end], - data=data['filteredMean'], - showDataPoint=options.showFittedPoint, - color=options.filteredColor, - label='filtered ' + data['name']) + plotData( + time=time[start:end], + data=data["filteredMean"], + showDataPoint=options.showFittedPoint, + color=options.filteredColor, + label="filtered " + data["name"], + ) if options.showConfidenceInterval: - upper, lower = getInterval(data['filteredMean'], - list(map(abs, data['filteredVar'])), - p=options.confidence) - - plotInterval(time=time[start:end], - upper=upper, lower=lower, - intervalType=options.intervalType, - color=options.filteredColor) - - plt.legend(loc='best', shadow=True) + upper, lower = getInterval( + data["filteredMean"], + list(map(abs, data["filteredVar"])), + p=options.confidence, + ) + + plotInterval( + time=time[start:end], + upper=upper, + lower=lower, + intervalType=options.intervalType, + color=options.filteredColor, + ) + + plt.legend(loc="best", shadow=True) # plot predicted results if needed if options.plotPredictedData: @@ -368,23 +456,30 @@ def plotComponent(time, data, result, options): start = result.filteredSteps[0] end = result.filteredSteps[1] + 1 - plotData(time=time[start:end], - data=data['predictedMean'], - showDataPoint=options.showFittedPoint, - color=options.predictedColor, - label='one-step predict for ' + data['name']) + plotData( + time=time[start:end], + data=data["predictedMean"], + showDataPoint=options.showFittedPoint, + color=options.predictedColor, + label="one-step predict for " + data["name"], + ) if options.showConfidenceInterval: - upper, lower = getInterval(data['predictedMean'], - list(map(abs, data['predictedVar'])), - p=options.confidence) - - plotInterval(time=time[start:end], - upper=upper, lower=lower, - intervalType=options.intervalType, - color=options.predictedColor) - - plt.legend(loc='best', shadow=True) + upper, lower = getInterval( + data["predictedMean"], + list(map(abs, data["predictedVar"])), + p=options.confidence, + ) + + plotInterval( + time=time[start:end], + upper=upper, + lower=lower, + intervalType=options.intervalType, + color=options.predictedColor, + ) + + plt.legend(loc="best", shadow=True) # plot smoothed results if needed if options.plotSmoothedData: @@ -394,28 +489,34 @@ def plotComponent(time, data, result, options): start = result.smoothedSteps[0] end = result.smoothedSteps[1] + 1 - plotData(time=time[start:end], - data=data['smoothedMean'], - showDataPoint=options.showFittedPoint, - color=options.smoothedColor, - label='smoothed ' + data['name']) + plotData( + time=time[start:end], + data=data["smoothedMean"], + showDataPoint=options.showFittedPoint, + color=options.smoothedColor, + label="smoothed " + data["name"], + ) if options.showConfidenceInterval: - upper, lower = getInterval(data['smoothedMean'], - list(map(abs, data['smoothedVar'])), - p=options.confidence) + upper, lower = getInterval( + data["smoothedMean"], + list(map(abs, data["smoothedVar"])), + p=options.confidence, + ) - plotInterval(time=time[start:end], - upper=upper, lower=lower, - intervalType=options.intervalType, - color=options.smoothedColor) + plotInterval( + time=time[start:end], + upper=upper, + lower=lower, + intervalType=options.intervalType, + color=options.smoothedColor, + ) - plt.legend(loc='best', shadow=True) + plt.legend(loc="best", shadow=True) # =========================== plot for prediction ======================== -def plotPrediction(time, data, predictedTime, - predictedData, predictedVar, options): +def plotPrediction(time, data, predictedTime, predictedData, predictedVar, options): """ Function for ploting N-day ahead prediction @@ -428,30 +529,38 @@ def plotPrediction(time, data, predictedTime, options: options for the plot, for details please refer to @dlm """ # plot the original data - plotData(time, data, showDataPoint=options.showDataPoint, - color=options.dataColor, - label='time series') + plotData( + time, + data, + showDataPoint=options.showDataPoint, + color=options.dataColor, + label="time series", + ) # overlay the predicted data - plotData(predictedTime, predictedData, - showDataPoint=options.showFittedPoint, - color=options.predictedColor, - label='N-day prediction') + plotData( + predictedTime, + predictedData, + showDataPoint=options.showFittedPoint, + color=options.predictedColor, + label="N-day prediction", + ) # overlay confidence interval if options.showConfidenceInterval: - upper, lower = getInterval(predictedData, - predictedVar, - p=options.confidence) - - plotInterval(time=predictedTime, - upper=upper, lower=lower, - intervalType=options.intervalType, - color=options.predictedColor) - - + upper, lower = getInterval(predictedData, predictedVar, p=options.confidence) + + plotInterval( + time=predictedTime, + upper=upper, + lower=lower, + intervalType=options.intervalType, + color=options.predictedColor, + ) + + # =========================== basic functions ============================ -def plotData(time, data, showDataPoint=True, color='black', label='unknown'): +def plotData(time, data, showDataPoint=True, color="black", label="unknown"): """ The function to plot data points. @@ -464,19 +573,19 @@ def plotData(time, data, showDataPoint=True, color='black', label='unknown'): """ if time is None: if showDataPoint: - plt.plot(data, 'o', color=color) - plt.plot(data, '-', color=color, label=label) + plt.plot(data, "o", color=color) + plt.plot(data, "-", color=color, label=label) else: - plt.plot(data, '-', color=color, label=label) + plt.plot(data, "-", color=color, label=label) else: if showDataPoint: - plt.plot(time, data, 'o', color=color) - plt.plot(time, data, '-', color=color, label=label) + plt.plot(time, data, "o", color=color) + plt.plot(time, data, "-", color=color, label=label) else: - plt.plot(time, data, '-', color=color, label=label) + plt.plot(time, data, "-", color=color, label=label) -def plotInterval(time, upper, lower, intervalType, color='black'): +def plotInterval(time, upper, lower, intervalType, color="black"): """ The function to plot confidence interval. @@ -488,19 +597,17 @@ def plotInterval(time, upper, lower, intervalType, color='black'): """ ALPHA = 0.4 if time is None: - if intervalType == 'line': - plt.plot(upper, '--', color=color) - plt.plot(lower, '--', color=color) - elif intervalType == 'ribbon': - plt.fill_between(upper, lower, facecolor=color, - alpha=ALPHA) + if intervalType == "line": + plt.plot(upper, "--", color=color) + plt.plot(lower, "--", color=color) + elif intervalType == "ribbon": + plt.fill_between(upper, lower, facecolor=color, alpha=ALPHA) else: - if intervalType == 'line': - plt.plot(time, upper, '--', color=color) - plt.plot(time, lower, '--', color=color) - elif intervalType == 'ribbon': - plt.fill_between(time, upper, lower, - facecolor=color, alpha=ALPHA) + if intervalType == "line": + plt.plot(time, upper, "--", color=color) + plt.plot(time, lower, "--", color=color) + elif intervalType == "ribbon": + plt.fill_between(time, upper, lower, facecolor=color, alpha=ALPHA) def plotInitialize(): diff --git a/pydlm/plot/dlmPlotMod.py b/pydlm/plot/dlmPlotMod.py index 23190c41..e1731e31 100644 --- a/pydlm/plot/dlmPlotMod.py +++ b/pydlm/plot/dlmPlotMod.py @@ -1,10 +1,12 @@ from pydlm.core._dlm import _dlm + class dlmPlotModule(_dlm): - """ A dlm module containing all plot methods. This is an API layer for the + """A dlm module containing all plot methods. This is an API layer for the `dlm` class. All methods defined in this class are public and can be called directly from `dlm` object. """ + def __init__(self, data, **options): super(dlmPlotModule, self).__init__(data, **options) @@ -16,7 +18,6 @@ def __init__(self, data, **options): # when plot function is called. self.plotLibLoaded = False - def turnOn(self, switch): """ "turn on" Operation for the dlm plotting options. @@ -31,29 +32,30 @@ def turnOn(self, switch): figures\n 'fitted dots', 'fitted' to plot fitted results with dots """ - if switch in set(['filtered plot', 'filter', - 'filtered results', 'filtering']): + if switch in set(["filtered plot", "filter", "filtered results", "filtering"]): self.options.plotFilteredData = True - elif switch in set(['smoothed plot', 'smooth', - 'smoothed results', 'smoothing']): + elif switch in set( + ["smoothed plot", "smooth", "smoothed results", "smoothing"] + ): self.options.plotSmoothedData = True - elif switch in set(['predict plot', 'predict', - 'predicted results', 'prediction']): + elif switch in set( + ["predict plot", "predict", "predicted results", "prediction"] + ): self.options.plotPredictedData = True - elif switch in set(['confidence interval', 'confidence', - 'interval', 'CI', 'ci']): + elif switch in set( + ["confidence interval", "confidence", "interval", "CI", "ci"] + ): self.options.showConfidenceInterval = True - elif switch in set(['data points', 'data point', 'points', 'data']): + elif switch in set(["data points", "data point", "points", "data"]): self.options.showDataPoint = True - elif switch in set(['multiple', 'multiple plots', - 'separate plots', 'separate']): + elif switch in set( + ["multiple", "multiple plots", "separate plots", "separate"] + ): self.options.separatePlot = True - elif switch in set(['fitted dots', 'fitted results', - 'fitted data', 'fitted']): + elif switch in set(["fitted dots", "fitted results", "fitted data", "fitted"]): self.options.showFittedPoint = True else: - raise NameError('no such options') - + raise NameError("no such options") def turnOff(self, switch): """ "turn off" Operation for the dlm plotting options. @@ -74,29 +76,30 @@ def turnOff(self, switch): results with dots """ - if switch in set(['filtered plot', 'filter', 'filtered results', - 'filtering']): + if switch in set(["filtered plot", "filter", "filtered results", "filtering"]): self.options.plotFilteredData = False - elif switch in set(['smoothed plot', 'smooth', 'smoothed results', - 'smoothing']): + elif switch in set( + ["smoothed plot", "smooth", "smoothed results", "smoothing"] + ): self.options.plotSmoothedData = False - elif switch in set(['predict plot', 'predict', 'predicted results', - 'prediction']): + elif switch in set( + ["predict plot", "predict", "predicted results", "prediction"] + ): self.options.plotPredictedData = False - elif switch in set(['confidence interval', 'confidence', 'interval', - 'CI', 'ci']): + elif switch in set( + ["confidence interval", "confidence", "interval", "CI", "ci"] + ): self.options.showConfidenceInterval = False - elif switch in set(['data points', 'data point', 'points', 'data']): + elif switch in set(["data points", "data point", "points", "data"]): self.options.showDataPoint = False - elif switch in set(['multiple', 'multiple plots', 'separate plots', - 'separate']): + elif switch in set( + ["multiple", "multiple plots", "separate plots", "separate"] + ): self.options.separatePlot = False - elif switch in set(['fitted dots', 'fitted results', 'fitted data', - 'fitted']): + elif switch in set(["fitted dots", "fitted results", "fitted data", "fitted"]): self.options.showFittedPoint = False else: - raise NameError('no such options') - + raise NameError("no such options") def setColor(self, switch, color): """ "set" Operation for the dlm plotting colors @@ -106,43 +109,35 @@ def setColor(self, switch, color): filtered/smoothed/predicted results, color: the color for the corresponding keyword. """ - if switch in set(['filtered plot', 'filter', 'filtered results', - 'filtering']): + if switch in set(["filtered plot", "filter", "filtered results", "filtering"]): self.options.filteredColor = color - elif switch in set(['smoothed plot', 'smooth', 'smoothed results', - 'smoothing']): + elif switch in set( + ["smoothed plot", "smooth", "smoothed results", "smoothing"] + ): self.options.smoothedColor = color - elif switch in set(['predict plot', 'predict', 'predicted results', - 'prediction']): + elif switch in set( + ["predict plot", "predict", "predicted results", "prediction"] + ): self.options.predictedColor = color - elif switch in set(['data points', 'data point', 'points', 'data']): + elif switch in set(["data points", "data point", "points", "data"]): self.options.dataColor = color else: - raise NameError('no such options') - + raise NameError("no such options") def setConfidence(self, p=0.95): - """ Set the confidence interval for the plot - - """ + """Set the confidence interval for the plot""" assert p >= 0 and p <= 1 self.options.confidence = p - def setIntervalType(self, intervalType): - """ Set the confidence interval type - - """ - if intervalType == 'ribbon' or intervalType == 'line': + """Set the confidence interval type""" + if intervalType == "ribbon" or intervalType == "line": self.options.intervalType = intervalType else: - raise NameError('No such type for confidence interval.') - + raise NameError("No such type for confidence interval.") def resetPlotOptions(self): - """ Reset the plotting option for the dlm class - - """ + """Reset the plotting option for the dlm class""" self.options.plotOriginalData = True self.options.plotFilteredData = True self.options.plotSmoothedData = True @@ -150,18 +145,17 @@ def resetPlotOptions(self): self.options.showDataPoint = False self.options.showFittedPoint = False self.options.showConfidenceInterval = True - self.options.dataColor = 'black' - self.options.filteredColor = 'blue' - self.options.predictedColor = 'green' - self.options.smoothedColor = 'red' + self.options.dataColor = "black" + self.options.filteredColor = "blue" + self.options.predictedColor = "green" + self.options.smoothedColor = "red" self.options.separatePlot = True self.options.confidence = 0.95 - self.options.intervalType = 'ribbon' - + self.options.intervalType = "ribbon" # plot the result according to the options - def plot(self, name='main'): - """ The main plot function. The dlmPlot and the matplotlib will only be loaded + def plot(self, name="main"): + """The main plot function. The dlmPlot and the matplotlib will only be loaded when necessary. Args: @@ -187,62 +181,60 @@ def plot(self, name='main'): # change option setting if some results are not available if not self.initialized: - raise NameError('The model must be constructed and' + - ' fitted before ploting.') + raise NameError( + "The model must be constructed and" + " fitted before ploting." + ) # check the filter status and automatically turn off bad plots self._checkPlotOptions() # plot the main time series after filtering - if name == 'main': + if name == "main": # if we just need one plot if self.options.separatePlot is not True: - dlmPlot.plotInOneFigure(time=time, - data=self.data, - result=self.result, - options=self.options) + dlmPlot.plotInOneFigure( + time=time, data=self.data, result=self.result, options=self.options + ) # otherwise, we plot in multiple figures else: - dlmPlot.plotInMultipleFigure(time=time, - data=self.data, - result=self.result, - options=self.options) + dlmPlot.plotInMultipleFigure( + time=time, data=self.data, result=self.result, options=self.options + ) # plot the component after filtering elif self._checkComponent(name): # create the data for ploting data = {} if self.options.plotFilteredData: - data['filteredMean'] = self.getMean( - filterType='forwardFilter', name=name) - data['filteredVar'] = self.getVar( - filterType='forwardFilter', name=name) + data["filteredMean"] = self.getMean( + filterType="forwardFilter", name=name + ) + data["filteredVar"] = self.getVar(filterType="forwardFilter", name=name) if self.options.plotSmoothedData: - data['smoothedMean'] = self.getMean( - filterType='backwardSmoother', name=name) - data['smoothedVar'] = self.getVar( - filterType='backwardSmoother', name=name) + data["smoothedMean"] = self.getMean( + filterType="backwardSmoother", name=name + ) + data["smoothedVar"] = self.getVar( + filterType="backwardSmoother", name=name + ) if self.options.plotPredictedData: - data['predictedMean'] = self.getMean( - filterType='predict', name=name) - data['predictedVar'] = self.getVar( - filterType='predict', name=name) + data["predictedMean"] = self.getMean(filterType="predict", name=name) + data["predictedVar"] = self.getVar(filterType="predict", name=name) if len(data) == 0: - raise NameError('Nothing is going to be drawn, due to ' + - 'user choices.') - data['name'] = name - dlmPlot.plotComponent(time=time, - data=data, - result=self.result, - options=self.options) + raise NameError( + "Nothing is going to be drawn, due to " + "user choices." + ) + data["name"] = name + dlmPlot.plotComponent( + time=time, data=data, result=self.result, options=self.options + ) dlmPlot.plotout() - def plotCoef(self, name, dimensions=None): - """ Plot function for the latent states (coefficents of dynamic + """Plot function for the latent states (coefficents of dynamic component). Args: @@ -267,15 +259,15 @@ def plotCoef(self, name, dimensions=None): # change option setting if some results are not available if not self.initialized: - raise NameError('The model must be constructed and' + - ' fitted before ploting.') + raise NameError( + "The model must be constructed and" + " fitted before ploting." + ) # check the filter status and automatically turn off bad plots self._checkPlotOptions() # plot the latent states for a given component if self._checkComponent(name): - # find its coordinates in the latent state indx = self.builder.componentIndex[name] # get the real latent states @@ -289,17 +281,18 @@ def plotCoef(self, name, dimensions=None): elif len(coordinates) > 5: coordinates = coordinates[:5] - dlmPlot.plotLatentState(time=time, - coordinates=coordinates, - result=self.result, - options=self.options, - name=name) + dlmPlot.plotLatentState( + time=time, + coordinates=coordinates, + result=self.result, + options=self.options, + name=name, + ) else: - raise NameError('No such component.') + raise NameError("No such component.") dlmPlot.plotout() - def plotPredictN(self, N=1, date=None, featureDict=None): """ Function to plot the N-day prediction results. @@ -333,28 +326,31 @@ def plotPredictN(self, N=1, date=None, featureDict=None): # change option setting if some results are not available if not self.initialized: - raise NameError('The model must be constructed and' + - ' fitted before ploting.') + raise NameError( + "The model must be constructed and" + " fitted before ploting." + ) # check the filter status and automatically turn off bad plots self._checkPlotOptions() predictedTimeRange = range(date, date + N) predictedData, predictedVar = self.predictN( - N=N, date=date, featureDict=featureDict) + N=N, date=date, featureDict=featureDict + ) dlmPlot.plotPrediction( - time=time, data=self.data, + time=time, + data=self.data, predictedTime=predictedTimeRange, predictedData=predictedData, predictedVar=predictedVar, - options=self.options) + options=self.options, + ) dlmPlot.plotout() - - def loadPlotLibrary(self): if not self.plotLibLoaded: global dlmPlot import pydlm.plot.dlmPlot as dlmPlot + self.plotLibLoaded = True diff --git a/pydlm/predict/_dlmPredict.py b/pydlm/predict/_dlmPredict.py index d0ac5def..959c5bb7 100644 --- a/pydlm/predict/_dlmPredict.py +++ b/pydlm/predict/_dlmPredict.py @@ -11,18 +11,17 @@ class _dlmPredict(_dlm): - """ The main class containing all prediction methods. + """The main class containing all prediction methods. Methods: _oneDayAheadPredict: predict one day a head. _continuePredict: continue predicting one day after _oneDayAheadPredict """ - # Note the following functions will modify the status of the model, so they # shall not be directly call through the main model if this behavior is not # desired. - + # featureDict contains all the features for prediction. # It is a dictionary with key equals to the name of the component and # the value as the new feature (a list). The function @@ -33,10 +32,10 @@ class _dlmPredict(_dlm): # (start_date, next_pred_date, [all_predicted_values]), which will be # used by _continuePredict. def _oneDayAheadPredict(self, date, featureDict=None): - """ One day ahead prediction based on the date and the featureDict. - The prediction could be on the last day and into the future or in + """One day ahead prediction based on the date and the featureDict. + The prediction could be on the last day and into the future or in the middle of the time series and ignore the rest. For predicting into - the future, the new features must be supplied to featureDict. For + the future, the new features must be supplied to featureDict. For prediction in the middle, the user can still supply the features which will be used priorily. The old features will be used if featureDict is None. @@ -56,15 +55,16 @@ def _oneDayAheadPredict(self, date, featureDict=None): A tuple of (predicted_mean, predicted_variance) """ if date > self.n - 1: - raise NameError('The date is beyond the data range.') + raise NameError("The date is beyond the data range.") # get the correct status of the model self._setModelStatus(date=date) self._constructEvaluationForPrediction( date=date + 1, featureDict=featureDict, - padded_data=self.padded_data[:(date + 1)]) - + padded_data=self.padded_data[: (date + 1)], + ) + # initialize the prediction status self.builder.model.prediction.step = 0 @@ -74,16 +74,15 @@ def _oneDayAheadPredict(self, date, featureDict=None): predictedObs = self.builder.model.prediction.obs predictedObsVar = self.builder.model.prediction.obsVar self.result.predictStatus = [ - date, # start_date - date + 1, # current_date - [predictedObs[0, 0]] # all historical predictions + date, # start_date + date + 1, # current_date + [predictedObs[0, 0]], # all historical predictions ] return (predictedObs, predictedObsVar) - def _continuePredict(self, featureDict=None): - """ Continue predicting one day after _oneDayAheadPredict or + """Continue predicting one day after _oneDayAheadPredict or after _continuePredict. After using _oneDayAheadPredict, the user can continue predicting by using _continuePredict. The featureDict act the same as in @@ -97,16 +96,18 @@ def _continuePredict(self, featureDict=None): A tuple of (predicted_mean, predicted_variance) """ if self.result.predictStatus is None: - raise NameError('_continoousPredict can only be used after ' + - '_oneDayAheadPredict') + raise NameError( + "_continoousPredict can only be used after " + "_oneDayAheadPredict" + ) startDate = self.result.predictStatus[0] currentDate = self.result.predictStatus[1] self._constructEvaluationForPrediction( date=currentDate + 1, featureDict=featureDict, - padded_data=self.padded_data[:(startDate + 1)] + - self.result.predictStatus[2]) + padded_data=self.padded_data[: (startDate + 1)] + + self.result.predictStatus[2], + ) self.Filter.predict(self.builder.model) predictedObs = self.builder.model.prediction.obs predictedObsVar = self.builder.model.prediction.obsVar @@ -114,13 +115,11 @@ def _continuePredict(self, featureDict=None): self.result.predictStatus[2].append(predictedObs[0, 0]) return (predictedObs, predictedObsVar) - # This function will modify the status of the object, use with caution. - def _constructEvaluationForPrediction(self, - date, - featureDict=None, - padded_data=None): - """ Construct the evaluation matrix based on date and featureDict. + def _constructEvaluationForPrediction( + self, date, featureDict=None, padded_data=None + ): + """Construct the evaluation matrix based on date and featureDict. Used for prediction. Features provided in the featureDict will be used preferrably. If the feature is not found in featureDict, the algorithm @@ -142,7 +141,7 @@ def _constructEvaluationForPrediction(self, if featureDict is not None: for name in featureDict: if name in self.builder.dynamicComponents: - comp = self.builder.dynamicComponents[name] + comp = self.builder.dynamicComponents[name] # the date is within range if date < comp.n: comp.features[date] = featureDict[name] @@ -151,7 +150,9 @@ def _constructEvaluationForPrediction(self, comp.features.append(featureDict[name]) comp.n += 1 else: - raise NameError("Feature is missing between the last predicted " + - "day and the new day") - + raise NameError( + "Feature is missing between the last predicted " + + "day and the new day" + ) + self.builder.updateEvaluation(date, padded_data) diff --git a/pydlm/predict/dlmPredictMod.py b/pydlm/predict/dlmPredictMod.py index ea2bfa7b..a81521e6 100644 --- a/pydlm/predict/dlmPredictMod.py +++ b/pydlm/predict/dlmPredictMod.py @@ -2,23 +2,23 @@ from numpy import matrix from pydlm.predict._dlmPredict import _dlmPredict + class dlmPredictModule(_dlmPredict): - """ A dlm module containing all prediction methods. This is an API layer + """A dlm module containing all prediction methods. This is an API layer for the `dlm` class. All methods defined in this class are public and can be called directly from `dlm` object. """ - # One day ahead prediction function def predict(self, date=None, featureDict=None): - """ One day ahead predict based on the current data. + """One day ahead predict based on the current data. The predict result is based on all the data before date and predict the - observation at date + days. + observation at date + days. - The prediction could be on the last day and into the future or in + The prediction could be on the last day and into the future or in the middle of the time series and ignore the rest. For predicting into - the future, the new features must be supplied to featureDict. For + the future, the new features must be supplied to featureDict. For prediction in the middle, the user can still supply the features which will be used priorily. The old features will be used if featureDict is None. @@ -43,19 +43,20 @@ def predict(self, date=None, featureDict=None): # check if the data on the date has been filtered if date > self.result.filteredSteps[1]: - raise NameError('Prediction can only be made right' + - ' after the filtered date') + raise NameError( + "Prediction can only be made right" + " after the filtered date" + ) # Clear the existing predictModel before the deepcopy to avoid recurrent # recurrent copy which could explode the memory and complexity. self._predictModel = None self._predictModel = deepcopy(self) - return self._predictModel._oneDayAheadPredict(date=date, - featureDict=featureDict) - + return self._predictModel._oneDayAheadPredict( + date=date, featureDict=featureDict + ) def continuePredict(self, featureDict=None): - """ Continue prediction after the one-day ahead predict. + """Continue prediction after the one-day ahead predict. If users want to have a multiple day prediction, they can opt to use continuePredict after predict with new features contained in @@ -79,17 +80,16 @@ def continuePredict(self, featureDict=None): A tupe. (predicted observation, variance) """ if self._predictModel is None: - raise NameError('continuePredict has to come after predict.') + raise NameError("continuePredict has to come after predict.") return self._predictModel._continuePredict(featureDict=featureDict) - # N day ahead prediction def predictN(self, N=1, date=None, featureDict=None): - """ N day ahead prediction based on the current data. + """N day ahead prediction based on the current data. This function is a convenient wrapper of predict() and - continuePredict(). If the prediction is into the future, i.e, > n, + continuePredict(). If the prediction is into the future, i.e, > n, the featureDict has to contain all feature vectors for multiple days for each dynamic component. For example, assume myDLM has a component named 'spy' which posseses two dimensions, @@ -114,7 +114,7 @@ def predictN(self, N=1, date=None, featureDict=None): """ if N < 1: - raise NameError('N has to be greater or equal to 1') + raise NameError("N has to be greater or equal to 1") # Take care if features are numpy matrix if featureDict is not None: for name in featureDict: @@ -124,8 +124,9 @@ def predictN(self, N=1, date=None, featureDict=None): predictedVar = [] # call predict for the first day - getSingleDayFeature = lambda f, i: ({k: v[i] for k, v in f.items()} - if f is not None else None) + getSingleDayFeature = lambda f, i: ( + {k: v[i] for k, v in f.items()} if f is not None else None + ) # Construct the single day featureDict featureDictOneDay = getSingleDayFeature(featureDict, 0) (obs, var) = self.predict(date=date, featureDict=featureDictOneDay) @@ -138,5 +139,7 @@ def predictN(self, N=1, date=None, featureDict=None): (obs, var) = self.continuePredict(featureDict=featureDictOneDay) predictedObs.append(obs) predictedVar.append(var) - return (self._1DmatrixToArray(predictedObs), - self._1DmatrixToArray(predictedVar)) + return ( + self._1DmatrixToArray(predictedObs), + self._1DmatrixToArray(predictedVar), + ) diff --git a/pydlm/tuner/_dlmTune.py b/pydlm/tuner/_dlmTune.py index 0492bd35..bb61b842 100644 --- a/pydlm/tuner/_dlmTune.py +++ b/pydlm/tuner/_dlmTune.py @@ -11,7 +11,7 @@ class _dlmTune(_dlm): - """ The main class containing all tuning methods. + """The main class containing all tuning methods. Methods: _getMSE: obtain the fitting model one-day ahead prediction MSE. @@ -19,33 +19,28 @@ class _dlmTune(_dlm): _setDiscounts: set discounts for different components. """ - # get the mse from the model def _getMSE(self): - if not self.initialized: - raise NameError('need to fit the model first') + raise NameError("need to fit the model first") if self.result.filteredSteps[1] == -1: - raise NameError('need to run forward filter first') + raise NameError("need to run forward filter first") mse = 0 - for i in range(self.result.filteredSteps[0], - self.result.filteredSteps[1] + 1): + for i in range(self.result.filteredSteps[0], self.result.filteredSteps[1] + 1): if self.data[i] is not None: mse += (self.data[i] - self.result.predictedObs[i]) ** 2 - mse = mse / (self.result.filteredSteps[1] + 1 - - self.result.filteredSteps[0]) - return mse[0,0] - + mse = mse / (self.result.filteredSteps[1] + 1 - self.result.filteredSteps[0]) + return mse[0, 0] # get the discount from the model def _getDiscounts(self): - if not self.initialized: - raise NameError('need to fit the model before one can' - 'fetch the discount factors') + raise NameError( + "need to fit the model before one can" "fetch the discount factors" + ) discounts = [] for comp in self.builder.componentIndex: @@ -53,23 +48,20 @@ def _getDiscounts(self): discounts.append(self.builder.discount[indx[0]]) return discounts - # set the model discount, this should never expose to the user # change the discount in the component would change the whole model. # change those in filter and builder only change the discounts # temporarily and will be corrected if we rebuild the model. def _setDiscounts(self, discounts, change_component=False): - if not self.initialized: - raise NameError('need to fit the model first') + raise NameError("need to fit the model first") for i, comp in enumerate(self.builder.componentIndex): indx = self.builder.componentIndex[comp] - self.builder.discount[indx[0]: (indx[1] + 1)] = discounts[i] + self.builder.discount[indx[0] : (indx[1] + 1)] = discounts[i] if change_component: component = self._fetchComponent(name=comp) - component.discount = self.builder.discount[indx[0]: (indx[1] + 1)] + component.discount = self.builder.discount[indx[0] : (indx[1] + 1)] self.Filter.updateDiscount(self.builder.discount) self.result.filteredSteps = [0, -1] - diff --git a/pydlm/tuner/dlmTuneMod.py b/pydlm/tuner/dlmTuneMod.py index d5df7250..d99ede58 100644 --- a/pydlm/tuner/dlmTuneMod.py +++ b/pydlm/tuner/dlmTuneMod.py @@ -3,13 +3,12 @@ import logging -class dlmTuneModule(_dlmTune): - """ A dlm model containing all tuning methods - """ +class dlmTuneModule(_dlmTune): + """A dlm model containing all tuning methods""" def getMSE(self): - """ Get the one-day ahead prediction mean square error. The mse is + """Get the one-day ahead prediction mean square error. The mse is estimated only for days that has been predicted. Returns: @@ -18,14 +17,13 @@ def getMSE(self): return self._getMSE() - def tune(self, maxit=100): - """ Automatic tuning of the discounting factors. + """Automatic tuning of the discounting factors. The method will call the model tuner class to use the default parameters to tune the discounting factors and change the discount factor permenantly. User needs to refit the model after tuning. - + If user wants a more refined tuning and not change any property of the existing model, they should opt to use the @modelTuner class. """ @@ -34,10 +32,10 @@ def tune(self, maxit=100): if self._logger.isEnabledFor(logging.INFO): self.fitForwardFilter() self._logger.info(f"The current mse is { str(self.getMSE()) }.") - + simpleTuner.tune(untunedDLM=self, maxit=maxit) self._setDiscounts(simpleTuner.getDiscounts(), change_component=True) - + if self._logger.isEnabledFor(logging.INFO): self.fitForwardFilter() self._logger.info(f"The new mse is { str(self.getMSE()) }.") diff --git a/pydlm/tuner/dlmTuner.py b/pydlm/tuner/dlmTuner.py index 953ba26e..791283f0 100644 --- a/pydlm/tuner/dlmTuner.py +++ b/pydlm/tuner/dlmTuner.py @@ -19,35 +19,36 @@ >>> mydlm.tune(maxit=100) This will permenantly change the discouting factor in mydlm. So if the user -prefer to build a new dlm with the new discount factor without changing the +prefer to build a new dlm with the new discount factor without changing the original one, one should opt to use the modelTuner class. """ + from copy import deepcopy from numpy import array import logging + class modelTuner: - """ The main class for modelTuner + """The main class for modelTuner Attributes: method: the optimization method. Currently only 'gradient_descent' is supported. loss: the optimization loss function. Currently only 'mse' (one-day ahead prediction) is supported. - - """ - def __init__(self, method='gradient_descent', loss='mse'): + """ + def __init__(self, method="gradient_descent", loss="mse"): self.method = method self.loss = loss self.current_mse = None self.err = 1e-4 self.discounts = None - def tune(self, untunedDLM, maxit=100, step = 1.0): - """ Main function for tuning the DLM model. + def tune(self, untunedDLM, maxit=100, step=1.0): + """Main function for tuning the DLM model. Args: untunedDLM: The DLM object that needs tuning @@ -64,13 +65,12 @@ def tune(self, untunedDLM, maxit=100, step = 1.0): tunedDLM.fitForwardFilter() discounts = array(tunedDLM._getDiscounts()) self.current_mse = tunedDLM._getMSE() - - # using gradient descent - if self.method == 'gradient_descent': + # using gradient descent + if self.method == "gradient_descent": # Disable all info and warning for faster processing. log_level = tunedDLM.getLoggingLevel() - tunedDLM.setLoggingLevel('CRITICAL') + tunedDLM.setLoggingLevel("CRITICAL") for i in range(maxit): gradient = self.find_gradient(discounts, tunedDLM) @@ -84,24 +84,31 @@ def tune(self, untunedDLM, maxit=100, step = 1.0): tunedDLM.setLoggingLevel(log_level) if i < maxit - 1: - tunedDLM._logger.info('Converge successfully!') + tunedDLM._logger.info("Converge successfully!") else: - tunedDLM._logger.warning('The algorithm stops without converging.') - if min(discounts) <= 0.7 + self.err or max(discounts) >= 1 - 2 * self.err: - tunedDLM._logger.info('Possible reason: some discount is too close to 1 or 0.7' - ' (0.7 is smallest discount that is permissible.') + tunedDLM._logger.warning("The algorithm stops without converging.") + if ( + min(discounts) <= 0.7 + self.err + or max(discounts) >= 1 - 2 * self.err + ): + tunedDLM._logger.info( + "Possible reason: some discount is too close to 1 or 0.7" + " (0.7 is smallest discount that is permissible." + ) else: - tunedDLM._logger.info('It might require more step to converge.' - ' Use tune(..., maixt = ) instead.') + tunedDLM._logger.info( + "It might require more step to converge." + " Use tune(..., maixt = ) instead." + ) self.discounts = discounts tunedDLM._setDiscounts(discounts, change_component=True) return tunedDLM def getDiscounts(self): - """ Get the tuned discounting factors. One for each component (even the - component being multi-dimensional, only one discounting factor will - be assigned to one component). Initialized to None. + """Get the tuned discounting factors. One for each component (even the + component being multi-dimensional, only one discounting factor will + be assigned to one component). Initialized to None. """ return self.discounts diff --git a/setup.py b/setup.py index 174a44f2..71bf7f5b 100644 --- a/setup.py +++ b/setup.py @@ -1,41 +1,43 @@ from setuptools import setup, find_packages from pathlib import Path + this_directory = Path(__file__).parent long_description = (this_directory / "README.md").read_text() setup( - name = 'pydlm', - author = 'Xiangyu Wang', - author_email = 'wwrechard@gmail.com', - description = ('A python library for the Bayesian dynamic linear ' + - 'model for time series modeling'), + name="pydlm", + author="Xiangyu Wang", + author_email="wwrechard@gmail.com", + description=( + "A python library for the Bayesian dynamic linear " + + "model for time series modeling" + ), long_description=long_description, - long_description_content_type='text/markdown', - license = 'BSD', - keywords = 'dlm bayes bayesian kalman filter smoothing dynamic model', - url = 'https://github.com/wwrechard/pydlm', - packages = find_packages(), - zip_safe= False, - classifiers = [ - 'Development Status :: 3 - Alpha', - 'Intended Audience :: Science/Research', - 'License :: OSI Approved :: BSD License', - 'Operating System :: OS Independent', - 'Programming Language :: Python', - 'Programming Language :: Python :: 3', - 'Topic :: Scientific/Engineering :: Artificial Intelligence' - ], - include_package_data = False, - install_requires = [ - 'numpy', - 'matplotlib', + long_description_content_type="text/markdown", + license="BSD", + keywords="dlm bayes bayesian kalman filter smoothing dynamic model", + url="https://github.com/wwrechard/pydlm", + packages=find_packages(), + zip_safe=False, + classifiers=[ + "Development Status :: 3 - Alpha", + "Intended Audience :: Science/Research", + "License :: OSI Approved :: BSD License", + "Operating System :: OS Independent", + "Programming Language :: Python", + "Programming Language :: Python :: 3", + "Topic :: Scientific/Engineering :: Artificial Intelligence", ], - tests_require = [ + include_package_data=False, + install_requires=[ + "numpy", + "matplotlib", ], - extras_require = { - 'docs': [ - 'Sphinx', + tests_require=[], + extras_require={ + "docs": [ + "Sphinx", ], }, ) diff --git a/tests/access/testDlmAccessMod.py b/tests/access/testDlmAccessMod.py index 4c5fae7e..a0b6e33e 100644 --- a/tests/access/testDlmAccessMod.py +++ b/tests/access/testDlmAccessMod.py @@ -10,161 +10,149 @@ class testDlmAccessMod(unittest.TestCase): def setUp(self): self.data5 = range(100) self.dlm5 = dlmAccessModule(self.data5) - self.dlm5.builder + trend(degree=0, discount=1, w=1.0) + \ - autoReg(degree=1, discount=1, w=1.0) + ( + self.dlm5.builder + + trend(degree=0, discount=1, w=1.0) + + autoReg(degree=1, discount=1, w=1.0) + ) self.dlm5._initialize() - self.dlm5.options.innovationType='whole' + self.dlm5.options.innovationType = "whole" def testGetMean(self): # for forward filter self.dlm5._forwardFilter(start=0, end=99, renew=False) self.dlm5.result.filteredSteps = [0, 99] - filteredTrend = self.dlm5.getMean(filterType='forwardFilter') + filteredTrend = self.dlm5.getMean(filterType="forwardFilter") self.assertEqual(len(filteredTrend), self.dlm5.n) np.testing.assert_array_equal( - filteredTrend, - self.dlm5._1DmatrixToArray(self.dlm5.result.filteredObs)) + filteredTrend, self.dlm5._1DmatrixToArray(self.dlm5.result.filteredObs) + ) # for predict filter - predictedTrend = self.dlm5.getMean(filterType='predict') + predictedTrend = self.dlm5.getMean(filterType="predict") np.testing.assert_array_equal( - predictedTrend, - self.dlm5._1DmatrixToArray(self.dlm5.result.predictedObs)) - + predictedTrend, self.dlm5._1DmatrixToArray(self.dlm5.result.predictedObs) + ) # for component with forward filter - arTrend = self.dlm5.getMean(filterType='forwardFilter', - name='ar2') + arTrend = self.dlm5.getMean(filterType="forwardFilter", name="ar2") trueAr = [item[1, 0] for item in self.dlm5.result.filteredState] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] - np.testing.assert_array_equal(arTrend, - self.dlm5._1DmatrixToArray(trueAr)) + np.testing.assert_array_equal(arTrend, self.dlm5._1DmatrixToArray(trueAr)) # for backward smoother self.dlm5._backwardSmoother(start=99) self.dlm5.result.smoothedSteps = [0, 99] - filteredTrend = self.dlm5.getMean(filterType='backwardSmoother') + filteredTrend = self.dlm5.getMean(filterType="backwardSmoother") np.testing.assert_array_equal( - filteredTrend, - self.dlm5._1DmatrixToArray(self.dlm5.result.smoothedObs)) + filteredTrend, self.dlm5._1DmatrixToArray(self.dlm5.result.smoothedObs) + ) # for component with backward smoother - arTrend = self.dlm5.getMean(filterType='backwardSmoother', - name='ar2') + arTrend = self.dlm5.getMean(filterType="backwardSmoother", name="ar2") trueAr = [item[1, 0] for item in self.dlm5.result.smoothedState] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] - np.testing.assert_array_equal(arTrend, - self.dlm5._1DmatrixToArray(trueAr)) - + np.testing.assert_array_equal(arTrend, self.dlm5._1DmatrixToArray(trueAr)) def testGetVar(self): # for forward filter self.dlm5._forwardFilter(start=0, end=99, renew=False) self.dlm5.result.filteredSteps = [0, 99] - filteredTrend = self.dlm5.getVar(filterType='forwardFilter') + filteredTrend = self.dlm5.getVar(filterType="forwardFilter") self.assertEqual(len(filteredTrend), self.dlm5.n) np.testing.assert_array_equal( - filteredTrend, - self.dlm5._1DmatrixToArray(self.dlm5.result.filteredObsVar)) + filteredTrend, self.dlm5._1DmatrixToArray(self.dlm5.result.filteredObsVar) + ) # for predict filter - predictedTrend = self.dlm5.getVar(filterType='predict') + predictedTrend = self.dlm5.getVar(filterType="predict") np.testing.assert_array_equal( - predictedTrend, - self.dlm5._1DmatrixToArray(self.dlm5.result.predictedObsVar)) + predictedTrend, self.dlm5._1DmatrixToArray(self.dlm5.result.predictedObsVar) + ) # for component with forward filter - arTrend = self.dlm5.getVar(filterType='forwardFilter', - name='ar2') + arTrend = self.dlm5.getVar(filterType="forwardFilter", name="ar2") trueAr = [item[1, 1] for item in self.dlm5.result.filteredCov] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] * comp.evaluation.T - np.testing.assert_array_equal(arTrend, - self.dlm5._1DmatrixToArray(trueAr)) + np.testing.assert_array_equal(arTrend, self.dlm5._1DmatrixToArray(trueAr)) # for backward smoother self.dlm5._backwardSmoother(start=99) self.dlm5.result.smoothedSteps = [0, 99] - filteredTrend = self.dlm5.getVar(filterType='backwardSmoother') + filteredTrend = self.dlm5.getVar(filterType="backwardSmoother") np.testing.assert_array_equal( - filteredTrend, - self.dlm5._1DmatrixToArray(self.dlm5.result.smoothedObsVar)) + filteredTrend, self.dlm5._1DmatrixToArray(self.dlm5.result.smoothedObsVar) + ) # for component with backward smoother - arTrend = self.dlm5.getVar(filterType='backwardSmoother', - name='ar2') + arTrend = self.dlm5.getVar(filterType="backwardSmoother", name="ar2") trueAr = [item[1, 1] for item in self.dlm5.result.smoothedCov] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] * comp.evaluation.T - np.testing.assert_array_equal(arTrend, - self.dlm5._1DmatrixToArray(trueAr)) - + np.testing.assert_array_equal(arTrend, self.dlm5._1DmatrixToArray(trueAr)) def testGetLatentState(self): # for forward filter self.dlm5._forwardFilter(start=0, end=99, renew=False) self.dlm5.result.filteredSteps = [0, 99] - filteredTrend = self.dlm5.getLatentState(filterType='forwardFilter') + filteredTrend = self.dlm5.getLatentState(filterType="forwardFilter") for i in range(100): np.testing.assert_array_equal( - filteredTrend[i], - self.dlm5.result.filteredState[i].flatten().tolist()) + filteredTrend[i], self.dlm5.result.filteredState[i].flatten().tolist() + ) # for predict filter - predictedTrend = self.dlm5.getLatentState(filterType='predict') + predictedTrend = self.dlm5.getLatentState(filterType="predict") for i in range(100): np.testing.assert_array_equal( - predictedTrend[i], - self.dlm5.result.predictedState[i].flatten().tolist()) - + predictedTrend[i], self.dlm5.result.predictedState[i].flatten().tolist() + ) + # for backward smoother self.dlm5._backwardSmoother(start=99) self.dlm5.result.smoothedSteps = [0, 99] - smoothedTrend = self.dlm5.getLatentState(filterType='backwardSmoother') + smoothedTrend = self.dlm5.getLatentState(filterType="backwardSmoother") for i in range(100): np.testing.assert_array_equal( - smoothedTrend[i], - self.dlm5.result.smoothedState[i].flatten().tolist()) - + smoothedTrend[i], self.dlm5.result.smoothedState[i].flatten().tolist() + ) def testGetLatentCov(self): # for forward filter self.dlm5._forwardFilter(start=0, end=99, renew=False) self.dlm5.result.filteredSteps = [0, 99] - filteredTrend = self.dlm5.getLatentCov(filterType='forwardFilter') - np.testing.assert_array_equal(filteredTrend, - self.dlm5.result.filteredCov) + filteredTrend = self.dlm5.getLatentCov(filterType="forwardFilter") + np.testing.assert_array_equal(filteredTrend, self.dlm5.result.filteredCov) # for predict filter - predictedTrend = self.dlm5.getLatentCov(filterType='predict') - np.testing.assert_array_equal(predictedTrend, - self.dlm5.result.predictedCov) + predictedTrend = self.dlm5.getLatentCov(filterType="predict") + np.testing.assert_array_equal(predictedTrend, self.dlm5.result.predictedCov) # for backward smoother self.dlm5._backwardSmoother(start=99) self.dlm5.result.smoothedSteps = [0, 99] - smoothedTrend = self.dlm5.getLatentCov(filterType='backwardSmoother') - np.testing.assert_array_equal(smoothedTrend, - self.dlm5.result.smoothedCov) + smoothedTrend = self.dlm5.getLatentCov(filterType="backwardSmoother") + np.testing.assert_array_equal(smoothedTrend, self.dlm5.result.smoothedCov) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/access/test_dlmGet.py b/tests/access/test_dlmGet.py index cc676740..5ef60a11 100644 --- a/tests/access/test_dlmGet.py +++ b/tests/access/test_dlmGet.py @@ -6,91 +6,92 @@ class test_dlmGet(unittest.TestCase): - - def setUp(self): self.data5 = range(100) self.dlm5 = _dlmGet(self.data5) - self.dlm5.builder + trend(degree=0, discount=1, w=1.0) + \ - autoReg(degree=1, discount=1, w=1.0) + ( + self.dlm5.builder + + trend(degree=0, discount=1, w=1.0) + + autoReg(degree=1, discount=1, w=1.0) + ) self.dlm5._initialize() - self.dlm5.options.innovationType='whole' - + self.dlm5.options.innovationType = "whole" def testGetLatentState(self): # for forward filter self.dlm5._forwardFilter(start=0, end=99, renew=False) self.dlm5.result.filteredSteps = [0, 99] filteredTrend = self.dlm5._getLatentState( - filterType='forwardFilter', name='trend', start=0, end=99) + filterType="forwardFilter", name="trend", start=0, end=99 + ) diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(filteredTrend[i][0] - - self.dlm5.result.filteredState[i][0, 0]) + diff += abs(filteredTrend[i][0] - self.dlm5.result.filteredState[i][0, 0]) self.assertAlmostEqual(diff, 0) # for prediction predictedTrend = self.dlm5._getLatentState( - filterType='predict', name='trend', start=0, end=99) + filterType="predict", name="trend", start=0, end=99 + ) diff = 0.0 for i in range(len(predictedTrend)): - diff += abs(predictedTrend[i][0] - - self.dlm5.result.predictedState[i][0, 0]) + diff += abs(predictedTrend[i][0] - self.dlm5.result.predictedState[i][0, 0]) self.assertAlmostEqual(diff, 0) # for backward smoother self.dlm5._backwardSmoother(start=99) self.dlm5.result.smoothedSteps = [0, 99] smoothedTrend = self.dlm5._getLatentState( - filterType='backwardSmoother', name='trend', start=0, end=99) + filterType="backwardSmoother", name="trend", start=0, end=99 + ) diff = 0.0 for i in range(len(smoothedTrend)): - diff += abs(smoothedTrend[i][0] - - self.dlm5.result.smoothedState[i][0, 0]) + diff += abs(smoothedTrend[i][0] - self.dlm5.result.smoothedState[i][0, 0]) self.assertAlmostEqual(diff, 0) - def testGetLatentCov(self): # for forward filter self.dlm5._forwardFilter(start=0, end=99, renew=False) self.dlm5.result.filteredSteps = [0, 99] filteredTrend = self.dlm5._getLatentCov( - filterType='forwardFilter', name='trend', start=0, end=99) + filterType="forwardFilter", name="trend", start=0, end=99 + ) diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(filteredTrend[i][0, 0] - - self.dlm5.result.filteredCov[i][0, 0]) + diff += abs(filteredTrend[i][0, 0] - self.dlm5.result.filteredCov[i][0, 0]) self.assertAlmostEqual(diff, 0) # for prediction predictedTrend = self.dlm5._getLatentCov( - filterType='predict', name='trend', start=0, end=99) + filterType="predict", name="trend", start=0, end=99 + ) diff = 0.0 for i in range(len(predictedTrend)): - diff += abs(predictedTrend[i][0, 0] - - self.dlm5.result.predictedCov[i][0, 0]) + diff += abs( + predictedTrend[i][0, 0] - self.dlm5.result.predictedCov[i][0, 0] + ) self.assertAlmostEqual(diff, 0) # for backward smoother self.dlm5._backwardSmoother(start=99) self.dlm5.result.smoothedSteps = [0, 99] smoothedTrend = self.dlm5._getLatentCov( - filterType='backwardSmoother', name='trend', start=0, end=99) + filterType="backwardSmoother", name="trend", start=0, end=99 + ) diff = 0.0 for i in range(len(smoothedTrend)): - diff += abs(smoothedTrend[i][0, 0] - - self.dlm5.result.smoothedCov[i][0, 0]) + diff += abs(smoothedTrend[i][0, 0] - self.dlm5.result.smoothedCov[i][0, 0]) self.assertAlmostEqual(diff, 0) - def testComponentMean(self): self.dlm5._forwardFilter(start=0, end=99, renew=False) self.dlm5.result.filteredSteps = [0, 99] # for component with forward filter - arTrend = self.dlm5._getComponentMean(filterType='forwardFilter', - name='ar2', start=0, end=99) + arTrend = self.dlm5._getComponentMean( + filterType="forwardFilter", name="ar2", start=0, end=99 + ) trueAr = [item[1, 0] for item in self.dlm5.result.filteredState] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] @@ -103,10 +104,11 @@ def testComponentMean(self): # for component with backward smoother self.dlm5._backwardSmoother(start=99) self.dlm5.result.smoothedSteps = [0, 99] - arTrend = self.dlm5._getComponentMean(filterType='backwardSmoother', - name='ar2', start=0, end=99) + arTrend = self.dlm5._getComponentMean( + filterType="backwardSmoother", name="ar2", start=0, end=99 + ) trueAr = [item[1, 0] for item in self.dlm5.result.smoothedState] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] @@ -116,15 +118,15 @@ def testComponentMean(self): diff += abs(arTrend[i] - trueAr[i]) self.assertAlmostEqual(diff, 0) - def testComponentVar(self): self.dlm5._forwardFilter(start=0, end=99, renew=False) self.dlm5.result.filteredSteps = [0, 99] # for component with forward filter - arTrend = self.dlm5._getComponentVar(filterType='forwardFilter', - name='ar2', start=0, end=99) + arTrend = self.dlm5._getComponentVar( + filterType="forwardFilter", name="ar2", start=0, end=99 + ) trueAr = [item[1, 1] for item in self.dlm5.result.filteredCov] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] * comp.evaluation.T @@ -137,10 +139,11 @@ def testComponentVar(self): # for component with backward smoother self.dlm5._backwardSmoother(start=99) self.dlm5.result.smoothedSteps = [0, 99] - arTrend = self.dlm5._getComponentVar(filterType='backwardSmoother', - name='ar2', start=0, end=99) + arTrend = self.dlm5._getComponentVar( + filterType="backwardSmoother", name="ar2", start=0, end=99 + ) trueAr = [item[1, 1] for item in self.dlm5.result.smoothedCov] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] * comp.evaluation.T @@ -151,5 +154,5 @@ def testComponentVar(self): self.assertAlmostEqual(diff, 0) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/base/testKalmanFilter.py b/tests/base/testKalmanFilter.py index 806c1409..373208f6 100644 --- a/tests/base/testKalmanFilter.py +++ b/tests/base/testKalmanFilter.py @@ -8,18 +8,15 @@ class testKalmanFilter(unittest.TestCase): - - def setUp(self): self.kf1 = kalmanFilter(discount=[1]) self.kf0 = kalmanFilter(discount=[1e-10]) self.kf11 = kalmanFilter(discount=[1, 1]) self.trend0 = trend(degree=0, discount=1, w=1.0) self.trend0_90 = trend(degree=0, discount=0.9, w=1.0) - self.trend0_98 = trend(degree=0, discount=0.98, w=1.0, name='a') + self.trend0_98 = trend(degree=0, discount=0.98, w=1.0, name="a") self.trend1 = trend(degree=1, discount=1, w=1.0) - def testForwardFilter(self): dlm = builder() dlm.add(self.trend0) @@ -53,7 +50,6 @@ def testForwardFilter(self): self.assertAlmostEqual(dlm.model.obs[0, 0], 1) self.assertAlmostEqual(dlm.model.prediction.obs[0, 0], 1) - def testForwardFilterMultiDim(self): dlm = builder() dlm.add(seasonality(period=2, discount=1, w=1.0)) @@ -67,7 +63,6 @@ def testForwardFilterMultiDim(self): self.assertAlmostEqual(dlm.model.state[0][0], -0.5) self.assertAlmostEqual(dlm.model.state[1][0], 0.5) - def testBackwardSmoother(self): dlm = builder() dlm.add(self.trend0) @@ -77,13 +72,10 @@ def testBackwardSmoother(self): # expect the smoothed mean at 1 will be 1/3, for discount = 1 self.kf1.forwardFilter(dlm.model, 1) self.kf1.forwardFilter(dlm.model, 0) - self.kf1.backwardSmoother(dlm.model, \ - np.array([[0.5]]), \ - np.array([[0.375]])) - self.assertAlmostEqual(dlm.model.obs[0, 0], 1.0/3) + self.kf1.backwardSmoother(dlm.model, np.array([[0.5]]), np.array([[0.375]])) + self.assertAlmostEqual(dlm.model.obs[0, 0], 1.0 / 3) self.assertAlmostEqual(dlm.model.sysVar[0, 0], 0.18518519) - # second order trend with discount = 1. The smoothed result should be # equal to a direct fit on the three data points, 0, 1, -1. Thus, the # smoothed observation should be 0.0 @@ -97,13 +89,10 @@ def testBackwardSmootherMultiDim(self): cov1 = dlm.model.sysVar self.kf11.forwardFilter(dlm.model, -1) - self.kf11.backwardSmoother(dlm.model, \ - rawState = state1, \ - rawSysVar = cov1) + self.kf11.backwardSmoother(dlm.model, rawState=state1, rawSysVar=cov1) self.assertAlmostEqual(dlm.model.obs[0, 0], 0.0) - def testMissingData(self): dlm = builder() dlm.add(self.trend0) @@ -115,40 +104,38 @@ def testMissingData(self): self.kf0.forwardFilter(dlm.model, None) self.assertAlmostEqual(dlm.model.obs[0, 0], 1.0) - self.assertAlmostEqual(dlm.model.obsVar[0, 0]/1e10, 0.5) + self.assertAlmostEqual(dlm.model.obsVar[0, 0] / 1e10, 0.5) self.kf0.forwardFilter(dlm.model, None) self.assertAlmostEqual(dlm.model.obs[0, 0], 1.0) - self.assertAlmostEqual(dlm.model.obsVar[0, 0]/1e10, 0.5) + self.assertAlmostEqual(dlm.model.obsVar[0, 0] / 1e10, 0.5) self.kf0.forwardFilter(dlm.model, 0) self.assertAlmostEqual(dlm.model.obs[0, 0], 0.0) - def testMissingEvaluation(self): dlm = builder() dlm.add(self.trend0) dlm.initialize() dlm.model.evaluation = np.array([[None]]) - self.kf1.forwardFilter(dlm.model, 1.0, dealWithMissingEvaluation = True) + self.kf1.forwardFilter(dlm.model, 1.0, dealWithMissingEvaluation=True) self.assertAlmostEqual(dlm.model.obs, 0.0) self.assertAlmostEqual(dlm.model.transition, 1.0) - def testEvolveMode(self): dlm = builder() dlm.add(self.trend0_90) dlm.add(self.trend0_98) dlm.initialize() - kf2 = kalmanFilter(discount=[0.9, 0.98], - updateInnovation='component', - index=dlm.componentIndex) + kf2 = kalmanFilter( + discount=[0.9, 0.98], updateInnovation="component", index=dlm.componentIndex + ) kf2.forwardFilter(dlm.model, 1.0) self.assertAlmostEqual(dlm.model.innovation[0, 1], 0.0) self.assertAlmostEqual(dlm.model.innovation[1, 0], 0.0) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/core/test_dlm.py b/tests/core/test_dlm.py index 833dc962..2ef01ed6 100644 --- a/tests/core/test_dlm.py +++ b/tests/core/test_dlm.py @@ -7,8 +7,6 @@ class test_dlm(unittest.TestCase): - - def setUp(self): self.data = [0] * 9 + [1] + [0] * 10 self.dlm1 = _dlm(self.data) @@ -17,87 +15,85 @@ def setUp(self): self.dlm2.builder + trend(degree=0, discount=1e-12, w=1.0) self.dlm1._initialize() self.dlm2._initialize() - self.dlm1.options.innovationType='whole' - self.dlm2.options.innovationType='whole' - + self.dlm1.options.innovationType = "whole" + self.dlm2.options.innovationType = "whole" def testForwardFilter(self): self.dlm1._forwardFilter(start=0, end=19, renew=False) self.assertAlmostEqual(np.sum(self.dlm1.result.filteredObs[0:9]), 0) - self.assertAlmostEqual(self.dlm1.result.filteredObs[9][0, 0], 1.0/11) - self.assertAlmostEqual(self.dlm1.result.filteredObs[19][0, 0], 1.0/21) + self.assertAlmostEqual(self.dlm1.result.filteredObs[9][0, 0], 1.0 / 11) + self.assertAlmostEqual(self.dlm1.result.filteredObs[19][0, 0], 1.0 / 21) self.dlm2._forwardFilter(start=0, end=19) self.assertAlmostEqual(np.sum(self.dlm2.result.filteredObs[0:9]), 0.0) self.assertAlmostEqual(self.dlm2.result.filteredObs[9][0, 0], 1.0) self.assertAlmostEqual(self.dlm2.result.filteredObs[19][0, 0], 0.0) - def testResetModelStatus(self): - self.dlm1._forwardFilter(start = 0, end = 19, renew = False) + self.dlm1._forwardFilter(start=0, end=19, renew=False) self.dlm1.result.filteredSteps = (0, 19) - self.assertAlmostEqual(self.dlm1.builder.model.obs, \ - self.dlm1.result.filteredObs[19]) + self.assertAlmostEqual( + self.dlm1.builder.model.obs, self.dlm1.result.filteredObs[19] + ) self.dlm1._resetModelStatus() - self.assertAlmostEqual(np.sum(self.dlm1.builder.model.state \ - - self.dlm1.builder.statePrior), 0.0) - + self.assertAlmostEqual( + np.sum(self.dlm1.builder.model.state - self.dlm1.builder.statePrior), 0.0 + ) def testSetModelStatus(self): - self.dlm1._forwardFilter(start = 0, end = 19, renew = False) + self.dlm1._forwardFilter(start=0, end=19, renew=False) self.dlm1.result.filteredSteps = (0, 19) - self.assertAlmostEqual(self.dlm1.builder.model.obs, \ - self.dlm1.result.filteredObs[19]) - self.dlm1._setModelStatus(date = 12) - self.assertAlmostEqual(self.dlm1.builder.model.obs, \ - self.dlm1.result.filteredObs[12]) - + self.assertAlmostEqual( + self.dlm1.builder.model.obs, self.dlm1.result.filteredObs[19] + ) + self.dlm1._setModelStatus(date=12) + self.assertAlmostEqual( + self.dlm1.builder.model.obs, self.dlm1.result.filteredObs[12] + ) def testForwaredFilterConsectiveness(self): - self.dlm1._forwardFilter(start = 0, end = 19, renew = False) + self.dlm1._forwardFilter(start=0, end=19, renew=False) filtered1 = self.dlm1.result.filteredObs self.dlm1._initialize() - self.dlm1._forwardFilter(start = 0, end = 13) + self.dlm1._forwardFilter(start=0, end=13) self.dlm1.result.filteredSteps = (0, 13) - self.dlm1._forwardFilter(start = 13, end = 19) + self.dlm1._forwardFilter(start=13, end=19) filtered2 = self.dlm1.result.filteredObs self.assertAlmostEqual(np.sum(np.array(filtered1) - np.array(filtered2)), 0.0) - def testBackwardSmoother(self): - self.dlm1._forwardFilter(start = 0, end = 19, renew = False) + self.dlm1._forwardFilter(start=0, end=19, renew=False) self.dlm1.result.filteredSteps = (0, 19) - self.dlm1._backwardSmoother(start = 19) - self.assertAlmostEqual(self.dlm1.result.smoothedObs[0][0, 0], 1.0/21) - self.assertAlmostEqual(self.dlm1.result.smoothedObs[19][0, 0], 1.0/21) + self.dlm1._backwardSmoother(start=19) + self.assertAlmostEqual(self.dlm1.result.smoothedObs[0][0, 0], 1.0 / 21) + self.assertAlmostEqual(self.dlm1.result.smoothedObs[19][0, 0], 1.0 / 21) - self.dlm2._forwardFilter(start = 0, end = 19) + self.dlm2._forwardFilter(start=0, end=19) self.dlm2.result.filteredSteps = (0, 19) - self.dlm2._backwardSmoother(start = 19) + self.dlm2._backwardSmoother(start=19) self.assertAlmostEqual(self.dlm2.result.smoothedObs[0][0, 0], 0.0) self.assertAlmostEqual(self.dlm2.result.smoothedObs[19][0, 0], 0.0) self.assertAlmostEqual(self.dlm2.result.smoothedObs[9][0, 0], 1.0) - def testLogger(self): - assert self.dlm1._logger == logging.getLogger('pydlm') - assert self.dlm2._logger == logging.getLogger('pydlm') + assert self.dlm1._logger == logging.getLogger("pydlm") + assert self.dlm2._logger == logging.getLogger("pydlm") assert self.dlm1._logger.getEffectiveLevel() == logging.INFO assert self.dlm2._logger.getEffectiveLevel() == logging.INFO def testSetAndGetLoggingLevel(self): - self.dlm1.setLoggingLevel('CRITICAL') + self.dlm1.setLoggingLevel("CRITICAL") assert self.dlm1._logger.getEffectiveLevel() == logging.CRITICAL - assert self.dlm1.getLoggingLevel() == 'CRITICAL' + assert self.dlm1.getLoggingLevel() == "CRITICAL" - self.dlm1.setLoggingLevel('INFO') + self.dlm1.setLoggingLevel("INFO") assert self.dlm1._logger.getEffectiveLevel() == logging.INFO - assert self.dlm1.getLoggingLevel() == 'INFO' + assert self.dlm1.getLoggingLevel() == "INFO" -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/modeler/testAutoReg.py b/tests/modeler/testAutoReg.py index e6143bb2..994b0bd2 100644 --- a/tests/modeler/testAutoReg.py +++ b/tests/modeler/testAutoReg.py @@ -3,21 +3,29 @@ class testAutoReg(unittest.TestCase): - - def setUp(self): self.data = [0, 1, 2, 3, 0, 1, 2, 3, 0, 1, 2, 3] - self.ar4 = autoReg(degree=4, name='ar4', padding=0, w=1.0) - + self.ar4 = autoReg(degree=4, name="ar4", padding=0, w=1.0) def testFeatureMatrix(self): - trueFeatures = [[0, 0, 0, 0], [0, 0, 0, 0], [0, 0, 0, 1], [0, 0, 1, 2], - [0, 1, 2, 3], [1, 2, 3, 0], [2, 3, 0, 1], [3, 0, 1, 2], - [0, 1, 2, 3], [1, 2, 3, 0], [2, 3, 0, 1], [3, 0, 1, 2]] + trueFeatures = [ + [0, 0, 0, 0], + [0, 0, 0, 0], + [0, 0, 0, 1], + [0, 0, 1, 2], + [0, 1, 2, 3], + [1, 2, 3, 0], + [2, 3, 0, 1], + [3, 0, 1, 2], + [0, 1, 2, 3], + [1, 2, 3, 0], + [2, 3, 0, 1], + [3, 0, 1, 2], + ] for i in range(12): self.ar4.updateEvaluation(i, self.data) self.assertEqual(self.ar4.evaluation.flatten().tolist(), trueFeatures[i]) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/modeler/testBuilder.py b/tests/modeler/testBuilder.py index 08764456..752d7949 100644 --- a/tests/modeler/testBuilder.py +++ b/tests/modeler/testBuilder.py @@ -9,26 +9,21 @@ class testBuilder(unittest.TestCase): - - def setUp(self): self.data = np.random.rand(10).tolist() self.features = np.random.rand(10, 2).tolist() self.trend = trend(degree=2, w=1.0) self.seasonality = seasonality(period=7, w=1.0) self.dynamic = dynamic(self.features, w=1.0) - self.autoReg = autoReg(degree=3, - w=1.0) + self.autoReg = autoReg(degree=3, w=1.0) self.builder1 = builder() self.builder2 = builder() - def testInitialization(self): self.assertEqual(len(self.builder1.staticComponents), 0) self.assertEqual(len(self.builder1.dynamicComponents), 0) self.assertEqual(len(self.builder1.automaticComponents), 0) - def testAddAndDelete(self): self.builder1 = self.builder1 + self.trend self.assertEqual(len(self.builder1.staticComponents), 1) @@ -45,69 +40,90 @@ def testAddAndDelete(self): self.assertEqual(len(self.builder1.dynamicComponents), 1) self.assertEqual(len(self.builder1.automaticComponents), 0) - self.builder1.delete('seasonality') + self.builder1.delete("seasonality") self.assertEqual(len(self.builder1.staticComponents), 1) self.assertEqual(len(self.builder1.dynamicComponents), 1) self.assertEqual(len(self.builder1.automaticComponents), 0) - self.assertEqual(self.builder1.staticComponents['trend'], self.trend) + self.assertEqual(self.builder1.staticComponents["trend"], self.trend) self.builder1 = self.builder1 + self.autoReg self.assertEqual(len(self.builder1.automaticComponents), 1) - def testInitialize(self): - self.builder1 = self.builder1 + self.trend + self.dynamic \ - + self.autoReg + self.builder1 = self.builder1 + self.trend + self.dynamic + self.autoReg self.builder1.initialize(data=self.data) - self.assertAlmostEqual(np.sum( - np.abs(self.builder1.model.evaluation - - mt.matrixAddByCol(mt.matrixAddByCol( - self.trend.evaluation, - self.dynamic.evaluation), - self.autoReg.evaluation))), 0.0) - + self.assertAlmostEqual( + np.sum( + np.abs( + self.builder1.model.evaluation + - mt.matrixAddByCol( + mt.matrixAddByCol( + self.trend.evaluation, self.dynamic.evaluation + ), + self.autoReg.evaluation, + ) + ) + ), + 0.0, + ) def testInitializeEvaluatoin(self): self.builder1 = self.builder1 + self.trend + self.dynamic - self.builder1.dynamicComponents['dynamic'].updateEvaluation(8) + self.builder1.dynamicComponents["dynamic"].updateEvaluation(8) self.builder1.initialize(data=self.data) - self.assertAlmostEqual(np.sum( - np.abs(self.builder1.model.evaluation - - mt.matrixAddByCol(self.trend.evaluation, - self.dynamic.evaluation))), 0.0) - + self.assertAlmostEqual( + np.sum( + np.abs( + self.builder1.model.evaluation + - mt.matrixAddByCol(self.trend.evaluation, self.dynamic.evaluation) + ) + ), + 0.0, + ) def testUpdate(self): - self.builder1 = self.builder1 + self.trend + self.dynamic \ - + self.autoReg + self.builder1 = self.builder1 + self.trend + self.dynamic + self.autoReg self.builder1.initialize(data=self.data) self.builder1.updateEvaluation(2, self.data) - self.assertAlmostEqual(np.sum( - np.abs(self.builder1.model.evaluation - - mt.matrixAddByCol(mt.matrixAddByCol( - self.trend.evaluation, - np.array([self.features[2]])), - self.autoReg.evaluation))), 0.0) - + self.assertAlmostEqual( + np.sum( + np.abs( + self.builder1.model.evaluation + - mt.matrixAddByCol( + mt.matrixAddByCol( + self.trend.evaluation, np.array([self.features[2]]) + ), + self.autoReg.evaluation, + ) + ) + ), + 0.0, + ) def testInitializFromBuilder(self): self.builder1 = self.builder1 + self.trend + self.dynamic - self.builder1.dynamicComponents['dynamic'].updateEvaluation(8) + self.builder1.dynamicComponents["dynamic"].updateEvaluation(8) self.builder1.initialize(data=self.data) self.data1 = self.data[8:-1] - self.builder2.initializeFromBuilder(exported_builder=self.builder1, data=self.data1) + self.builder2.initializeFromBuilder( + exported_builder=self.builder1, data=self.data1 + ) # make sure builder1 and builder2 are identical. - self.assertDictEqual(self.builder1.staticComponents, - self.builder2.staticComponents) - self.assertDictEqual(self.builder1.dynamicComponents, - self.builder2.dynamicComponents) - self.assertDictEqual(self.builder1.automaticComponents, - self.builder2.automaticComponents) - -if __name__ == '__main__': + self.assertDictEqual( + self.builder1.staticComponents, self.builder2.staticComponents + ) + self.assertDictEqual( + self.builder1.dynamicComponents, self.builder2.dynamicComponents + ) + self.assertDictEqual( + self.builder1.automaticComponents, self.builder2.automaticComponents + ) + + +if __name__ == "__main__": unittest.main() diff --git a/tests/modeler/testComponent.py b/tests/modeler/testComponent.py new file mode 100644 index 00000000..76ab6fa9 --- /dev/null +++ b/tests/modeler/testComponent.py @@ -0,0 +1,37 @@ +from pydlm.modeler.component import component + +import unittest + + +class testComponent(unittest.TestCase): + def testCreation(self): + test_component = component() + + self.assertTrue(hasattr(test_component, "d")) + self.assertTrue(hasattr(test_component, "name")) + self.assertTrue(hasattr(test_component, "componentType")) + self.assertTrue(hasattr(test_component, "discount")) + self.assertTrue(hasattr(test_component, "evaluation")) + self.assertTrue(hasattr(test_component, "transition")) + self.assertTrue(hasattr(test_component, "covPrior")) + self.assertTrue(hasattr(test_component, "meanPrior")) + + def testEqual(self): + component_1 = component() + component_2 = component() + + component_1.d = 1 + component_1.name = "component_1" + component_1.discount = 1 + + component_2.d = 1 + component_2.name = "component_2" + component_2.discount = 1 + self.assertTrue(component_1 != component_2) + + component_2.name = "component_1" + self.assertTrue(component_1 == component_2) + + +if __name__ == "__main__": + unittest.main() diff --git a/tests/modeler/testDynamic.py b/tests/modeler/testDynamic.py index d29c089e..2448142a 100644 --- a/tests/modeler/testDynamic.py +++ b/tests/modeler/testDynamic.py @@ -4,20 +4,16 @@ class testDynamic(unittest.TestCase): - - def setUp(self): self.features = np.random.rand(10, 2).tolist() self.features2 = np.random.rand(10, 1).tolist() self.newDynamic = dynamic(features=self.features, w=1.0) self.newDynamic2 = dynamic(features=self.features2, w=1.0) - def testInputNumpyMatrix(self): dynamic(features=np.random.rand(10, 2), w=1.0) pass - def testInitialization(self): self.assertEqual(self.newDynamic.d, 2) self.assertEqual(self.newDynamic.n, 10) @@ -27,34 +23,31 @@ def testInitialization(self): self.assertEqual(self.newDynamic2.n, 10) self.assertEqual(self.newDynamic2.features, self.features2) - - def testUpdate(self): for i in range(10): self.newDynamic.updateEvaluation(i) np.testing.assert_array_equal( - self.newDynamic.evaluation, np.array([self.features[i]])) + self.newDynamic.evaluation, np.array([self.features[i]]) + ) for i in range(10): self.newDynamic2.updateEvaluation(i) np.testing.assert_array_equal( - self.newDynamic2.evaluation, np.array([self.features2[i]])) - + self.newDynamic2.evaluation, np.array([self.features2[i]]) + ) def testAppendNewData(self): self.newDynamic.appendNewData([[1, 2]]) self.assertEqual(self.newDynamic.features[-1], [1, 2]) - def testPopout(self): self.newDynamic.popout(0) self.assertEqual(self.newDynamic.features, self.features[1:]) - def testAlter(self): self.newDynamic.alter(1, [0, 0]) self.assertEqual(self.newDynamic.features[1], [0, 0]) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/modeler/testLongSeason.py b/tests/modeler/testLongSeason.py index 41447b0e..f04146e8 100644 --- a/tests/modeler/testLongSeason.py +++ b/tests/modeler/testLongSeason.py @@ -3,45 +3,62 @@ class testLongSeason(unittest.TestCase): - - def setUp(self): self.longSeason = longSeason(period=4, stay=4, w=1.0) self.longSeason2 = longSeason(period=2, stay=3, w=1.0) - def testLongSeasonProperties(self): - self.assertEqual(self.longSeason.componentType, 'longSeason') + self.assertEqual(self.longSeason.componentType, "longSeason") self.assertEqual(self.longSeason.period, 4) self.assertEqual(self.longSeason.stay, 4) - self.assertEqual(self.longSeason2.componentType, 'longSeason') + self.assertEqual(self.longSeason2.componentType, "longSeason") self.assertEqual(self.longSeason2.period, 2) self.assertEqual(self.longSeason2.stay, 3) - def testUpdateEvaluation(self): - trueFeatures = [[1, 0, 0, 0], [1, 0, 0, 0], [1, 0, 0, 0], [1, 0, 0, 0], - [0, 1, 0, 0], [0, 1, 0, 0], [0, 1, 0, 0], [0, 1, 0, 0], - [0, 0, 1, 0], [0, 0, 1, 0], [0, 0, 1, 0], [0, 0, 1, 0], - [0, 0, 0, 1]] + trueFeatures = [ + [1, 0, 0, 0], + [1, 0, 0, 0], + [1, 0, 0, 0], + [1, 0, 0, 0], + [0, 1, 0, 0], + [0, 1, 0, 0], + [0, 1, 0, 0], + [0, 1, 0, 0], + [0, 0, 1, 0], + [0, 0, 1, 0], + [0, 0, 1, 0], + [0, 0, 1, 0], + [0, 0, 0, 1], + ] for i in range(13): self.longSeason.updateEvaluation(step=i) - self.assertEqual(self.longSeason.evaluation.flatten().tolist(), - trueFeatures[i]) + self.assertEqual( + self.longSeason.evaluation.flatten().tolist(), trueFeatures[i] + ) def testUpdateEvaluation2(self): - trueFeatures = [[1, 0], [1, 0], [1, 0], - [0, 1], [0, 1], [0, 1], - [1, 0], [1, 0], [1, 0], - [0, 1], [0, 1], [0, 1], - ] + trueFeatures = [ + [1, 0], + [1, 0], + [1, 0], + [0, 1], + [0, 1], + [0, 1], + [1, 0], + [1, 0], + [1, 0], + [0, 1], + [0, 1], + [0, 1], + ] for i in range(12): self.longSeason2.updateEvaluation(step=i) - self.assertEqual(self.longSeason2.evaluation.flatten().tolist(), - trueFeatures[i]) - + self.assertEqual( + self.longSeason2.evaluation.flatten().tolist(), trueFeatures[i] + ) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/modeler/testMatrixTools.py b/tests/modeler/testMatrixTools.py index 83c6c3df..f5851135 100644 --- a/tests/modeler/testMatrixTools.py +++ b/tests/modeler/testMatrixTools.py @@ -3,8 +3,8 @@ from pydlm.modeler.matrixTools import matrixTools as mt -class testMatrixTools(unittest.TestCase): +class testMatrixTools(unittest.TestCase): def testMatrixAddByRow(self): A = np.array([[1], [2]]) B = np.array([[3]]) @@ -41,6 +41,7 @@ def addTwoVectors(self): np.testing.assert_array_equal(mt.addTwoVectors(a, None), a) np.testing.assert_array_equal(mt.addTwoVectors(None, b), b) - -if __name__ == '__main__': + + +if __name__ == "__main__": unittest.main() diff --git a/tests/modeler/testSeasonality.py b/tests/modeler/testSeasonality.py index 2cc21d74..94cc19fc 100644 --- a/tests/modeler/testSeasonality.py +++ b/tests/modeler/testSeasonality.py @@ -3,19 +3,17 @@ from pydlm.modeler.seasonality import seasonality -class testSeasonality(unittest.TestCase): - +class testSeasonality(unittest.TestCase): def testInitialization(self): newSeasonality = seasonality(period=7) newSeasonality.checkDimensions() - def testFreeForm(self): seasonality2 = seasonality(period=2, discount=1, w=1.0) self.assertEqual(seasonality2.meanPrior.tolist(), [[0], [0]]) self.assertEqual(seasonality2.covPrior.tolist(), [[0.5, -0.5], [-0.5, 0.5]]) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/modeler/testTrends.py b/tests/modeler/testTrends.py index 0f5beaee..24d8948c 100644 --- a/tests/modeler/testTrends.py +++ b/tests/modeler/testTrends.py @@ -1,17 +1,15 @@ import pydlm import unittest -class testTrend(unittest.TestCase): - +class testTrend(unittest.TestCase): def setUp(self): self.DEGREE = 3 - def testInitialization(self): newTrend = pydlm.modeler.trends.trend(self.DEGREE) newTrend.checkDimensions() -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/plot/testDlmPlot.py b/tests/plot/testDlmPlot.py index a2c0fdc4..b35182e0 100644 --- a/tests/plot/testDlmPlot.py +++ b/tests/plot/testDlmPlot.py @@ -3,12 +3,13 @@ from pydlm.modeler.trends import trend from pydlm.dlm import dlm -class testDlmPlot(unittest.TestCase): +class testDlmPlot(unittest.TestCase): def testPlot(self): dlm1 = dlm(range(100)) + trend(1) dlm1.fit() dlm1.plot() -if __name__ == '__main__': + +if __name__ == "__main__": unittest.main() diff --git a/tests/predict/test_dlmPredict.py b/tests/predict/test_dlmPredict.py index 2a379d61..bcd98581 100644 --- a/tests/predict/test_dlmPredict.py +++ b/tests/predict/test_dlmPredict.py @@ -8,8 +8,6 @@ class test_dlmPredict(unittest.TestCase): - - def setUp(self): self.data = [0] * 9 + [1] + [0] * 10 self.data5 = range(100) @@ -18,70 +16,66 @@ def setUp(self): self.dlm5 = _dlmPredict(self.data5) self.dlm3.builder + seasonality(period=2, discount=1, w=1.0) - self.dlm4.builder + dynamic(features=[[0] for i in range(5)] + - [[1] for i in range(5)], discount=1, - w=1.0) - self.dlm5.builder + trend(degree=0, discount=1, w=1.0) + \ - autoReg(degree=1, discount=1, w=1.0) + self.dlm4.builder + dynamic( + features=[[0] for i in range(5)] + [[1] for i in range(5)], + discount=1, + w=1.0, + ) + ( + self.dlm5.builder + + trend(degree=0, discount=1, w=1.0) + + autoReg(degree=1, discount=1, w=1.0) + ) self.dlm3._initialize() self.dlm4._initialize() self.dlm5._initialize() - self.dlm3.options.innovationType='whole' - self.dlm4.options.innovationType='whole' - self.dlm5.options.innovationType='whole' - + self.dlm3.options.innovationType = "whole" + self.dlm4.options.innovationType = "whole" + self.dlm5.options.innovationType = "whole" def testOneDayAheadPredictWithoutDynamic(self): self.dlm3._forwardFilter(start=0, end=11, renew=False) self.dlm3.result.filteredSteps = (0, 11) (obs, var) = self.dlm3._oneDayAheadPredict(date=11) - self.assertAlmostEqual(obs, -6.0/7) - self.assertAlmostEqual(self.dlm3.result.predictStatus, - [11, 12, [-6.0/7]]) + self.assertAlmostEqual(obs, -6.0 / 7) + self.assertAlmostEqual(self.dlm3.result.predictStatus, [11, 12, [-6.0 / 7]]) (obs, var) = self.dlm3._oneDayAheadPredict(date=2) - self.assertAlmostEqual(obs, 3.0/5) + self.assertAlmostEqual(obs, 3.0 / 5) # notice that the two latent states always sum up to 0 - self.assertAlmostEqual(self.dlm3.result.predictStatus, - [2, 3, [3.0/5]]) - + self.assertAlmostEqual(self.dlm3.result.predictStatus, [2, 3, [3.0 / 5]]) def testOneDayAheadPredictWithDynamic(self): self.dlm4._forwardFilter(start=0, end=9, renew=False) self.dlm4.result.filteredSteps = (0, 9) - featureDict = {'dynamic': 2.0} - (obs, var) = self.dlm4._oneDayAheadPredict(date=9, - featureDict=featureDict) - self.assertAlmostEqual(obs, 5.0/6 * 2) - + featureDict = {"dynamic": 2.0} + (obs, var) = self.dlm4._oneDayAheadPredict(date=9, featureDict=featureDict) + self.assertAlmostEqual(obs, 5.0 / 6 * 2) def testContinuePredictWithoutDynamic(self): self.dlm3._forwardFilter(start=0, end=11, renew=False) self.dlm3.result.filteredSteps = (0, 11) (obs, var) = self.dlm3._oneDayAheadPredict(date=11) - self.assertAlmostEqual(self.dlm3.result.predictStatus, - [11, 12, [-6.0/7]]) + self.assertAlmostEqual(self.dlm3.result.predictStatus, [11, 12, [-6.0 / 7]]) (obs, var) = self.dlm3._continuePredict() - self.assertAlmostEqual(self.dlm3.result.predictStatus, - [11, 13, [-6.0/7, 6.0/7]]) - + self.assertAlmostEqual( + self.dlm3.result.predictStatus, [11, 13, [-6.0 / 7, 6.0 / 7]] + ) def testContinuePredictWithDynamic(self): self.dlm4._forwardFilter(start=0, end=9, renew=False) self.dlm4.result.filteredSteps = (0, 9) - featureDict = {'dynamic': 2.0} - (obs, var) = self.dlm4._oneDayAheadPredict(date=9, - featureDict=featureDict) - self.assertAlmostEqual(self.dlm4.result.predictStatus, - [9, 10, [5.0/6 * 2]]) + featureDict = {"dynamic": 2.0} + (obs, var) = self.dlm4._oneDayAheadPredict(date=9, featureDict=featureDict) + self.assertAlmostEqual(self.dlm4.result.predictStatus, [9, 10, [5.0 / 6 * 2]]) - featureDict = {'dynamic': 3.0} + featureDict = {"dynamic": 3.0} (obs, var) = self.dlm4._continuePredict(featureDict=featureDict) - self.assertAlmostEqual(self.dlm4.result.predictStatus, - [9, 11, [5.0/6 * 2, 5.0/6 * 3]]) - + self.assertAlmostEqual( + self.dlm4.result.predictStatus, [9, 11, [5.0 / 6 * 2, 5.0 / 6 * 3]] + ) def testPredictWithAutoReg(self): self.dlm5._forwardFilter(start=0, end=99, renew=False) @@ -92,5 +86,5 @@ def testPredictWithAutoReg(self): self.assertAlmostEqual(obs[0, 0], 101.07480945) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/testDlm.py b/tests/testDlm.py index 2b6027a6..ff93e55a 100644 --- a/tests/testDlm.py +++ b/tests/testDlm.py @@ -7,9 +7,8 @@ from pydlm.modeler.autoReg import autoReg from pydlm.dlm import dlm -class testDlm(unittest.TestCase): - +class testDlm(unittest.TestCase): def setUp(self): self.data = [0] * 9 + [1] + [0] * 10 self.data5 = range(100) @@ -25,60 +24,65 @@ def setUp(self): self.dlm1 + trend(degree=0, discount=1, w=1.0) self.dlm2 + trend(degree=0, discount=1e-12, w=1.0) self.dlm3 + seasonality(period=2, discount=1, w=1.0) - self.dlm4 + dynamic(features=[[0] for i in range(5)] + - [[1] for i in range(5)], discount=1, w=1.0) - self.dlm5 + trend(degree=0, discount=1, w=1.0) + \ - autoReg(degree=1, discount=1, w=1.0) - self.dlm6 + trend(degree=0, discount=1, w=1.0) + \ - autoReg(degree=2, discount=1, w=1.0) - self.dlm1.evolveMode('dependent') - self.dlm2.evolveMode('dependent') - self.dlm3.evolveMode('dependent') - self.dlm4.evolveMode('dependent') - self.dlm5.evolveMode('dependent') - self.dlm6.evolveMode('dependent') - + self.dlm4 + dynamic( + features=[[0] for i in range(5)] + [[1] for i in range(5)], + discount=1, + w=1.0, + ) + ( + self.dlm5 + + trend(degree=0, discount=1, w=1.0) + + autoReg(degree=1, discount=1, w=1.0) + ) + ( + self.dlm6 + + trend(degree=0, discount=1, w=1.0) + + autoReg(degree=2, discount=1, w=1.0) + ) + self.dlm1.evolveMode("dependent") + self.dlm2.evolveMode("dependent") + self.dlm3.evolveMode("dependent") + self.dlm4.evolveMode("dependent") + self.dlm5.evolveMode("dependent") + self.dlm6.evolveMode("dependent") def testAdd(self): - trend2 = trend(2, name='trend2') + trend2 = trend(2, name="trend2") self.dlm1 = self.dlm1 + trend2 - self.assertEqual(self.dlm1.builder.staticComponents['trend2'], trend2) + self.assertEqual(self.dlm1.builder.staticComponents["trend2"], trend2) - dynamic2 = dynamic(features=self.features, name='d2') + dynamic2 = dynamic(features=self.features, name="d2") self.dlm1 = self.dlm1 + dynamic2 - self.assertEqual(self.dlm1.builder.dynamicComponents['d2'], dynamic2) + self.assertEqual(self.dlm1.builder.dynamicComponents["d2"], dynamic2) - ar3 = autoReg(degree=3, name='ar3') + ar3 = autoReg(degree=3, name="ar3") self.dlm1 = self.dlm1 + ar3 - self.assertEqual(self.dlm1.builder.automaticComponents['ar3'], ar3) - + self.assertEqual(self.dlm1.builder.automaticComponents["ar3"], ar3) def testDelete(self): - trend2 = trend(2, name='trend2') + trend2 = trend(2, name="trend2") self.dlm1 = self.dlm1 + trend2 - self.dlm1.delete('trend2') + self.dlm1.delete("trend2") self.assertEqual(len(self.dlm1.builder.staticComponents), 1) - def testFitForwardFilter(self): - self.dlm1.fitForwardFilter(useRollingWindow = False) + self.dlm1.fitForwardFilter(useRollingWindow=False) self.assertEqual(self.dlm1.result.filteredSteps, [0, 19]) self.assertAlmostEqual(np.sum(self.dlm1.result.filteredObs[0:9]), 0) - self.assertAlmostEqual(self.dlm1.result.filteredObs[9][0, 0], 1.0/11) - self.assertAlmostEqual(self.dlm1.result.filteredObs[19][0, 0], 1.0/21) + self.assertAlmostEqual(self.dlm1.result.filteredObs[9][0, 0], 1.0 / 11) + self.assertAlmostEqual(self.dlm1.result.filteredObs[19][0, 0], 1.0 / 21) - self.dlm2.fitForwardFilter(useRollingWindow = False) + self.dlm2.fitForwardFilter(useRollingWindow=False) self.assertAlmostEqual(np.sum(self.dlm2.result.filteredObs[0:9]), 0.0) self.assertAlmostEqual(self.dlm2.result.filteredObs[9][0, 0], 1.0) self.assertAlmostEqual(self.dlm2.result.filteredObs[19][0, 0], 0.0) - def testFitBackwardSmoother(self): self.dlm1.fitForwardFilter() self.dlm1.fitBackwardSmoother() self.assertEqual(self.dlm1.result.smoothedSteps, [0, 19]) - self.assertAlmostEqual(self.dlm1.result.smoothedObs[0][0, 0], 1.0/21) - self.assertAlmostEqual(self.dlm1.result.smoothedObs[19][0, 0], 1.0/21) + self.assertAlmostEqual(self.dlm1.result.smoothedObs[0][0, 0], 1.0 / 21) + self.assertAlmostEqual(self.dlm1.result.smoothedObs[19][0, 0], 1.0 / 21) self.dlm2.fitForwardFilter() self.dlm2.fitBackwardSmoother() @@ -86,45 +90,61 @@ def testFitBackwardSmoother(self): self.assertAlmostEqual(self.dlm2.result.smoothedObs[19][0, 0], 0.0) self.assertAlmostEqual(self.dlm2.result.smoothedObs[9][0, 0], 1.0) - def testAppend(self): dlm4 = dlm(self.data[0:11]) dlm4 + self.trend0 - dlm4.evolveMode('dependent') + dlm4.evolveMode("dependent") dlm4.fitForwardFilter() self.assertEqual(dlm4.n, 11) - dlm4.append(self.data[11 : 20]) + dlm4.append(self.data[11:20]) self.assertEqual(dlm4.n, 20) dlm4.fitForwardFilter() self.dlm1.fitForwardFilter() - self.assertAlmostEqual(np.sum(np.array(dlm4.result.filteredObs) - - np.array(self.dlm1.result.filteredObs)), 0.0) - + self.assertAlmostEqual( + np.sum( + np.array(dlm4.result.filteredObs) + - np.array(self.dlm1.result.filteredObs) + ), + 0.0, + ) + + def testAppendEmptyDLM(self): + dlm_original = dlm(self.data[0:20]) + dlm_original + trend(degree=1, discount=1) + autoReg(degree=2, discount=1) + dlm_original.fitForwardFilter() + + dlm_append = dlm(data=[]) + dlm_append + trend(degree=1, discount=1) + autoReg(degree=2, discount=1) + dlm_append.append(self.data[0:20]) + dlm_append.fitForwardFilter() + + np.testing.assert_array_equal(dlm_original.getMean(), dlm_append.getMean()) def testAppendDynamic(self): # we feed the data to dlm4 via two segments dlm4 = dlm(self.data[0:11]) - dlm4 + self.trend1 + dynamic(features = self.features[0:11], - discount = 1) + dlm4 + self.trend1 + dynamic(features=self.features[0:11], discount=1) dlm4.fitForwardFilter() - dlm4.append(self.data[11 : 20]) - dlm4.append(self.features[11 : 20], component = 'dynamic') + dlm4.append(self.data[11:20]) + dlm4.append(self.features[11:20], component="dynamic") dlm4.fitForwardFilter() # we feed the data to dlm5 all at once dlm5 = dlm(self.data) - dlm5 + self.trend1 + dynamic(features = self.features, - discount = 1) + dlm5 + self.trend1 + dynamic(features=self.features, discount=1) dlm5.fitForwardFilter() - self.assertAlmostEqual(np.sum(np.array(dlm4.result.filteredObs) - - np.array(dlm5.result.filteredObs)), 0.0) - + self.assertAlmostEqual( + np.sum( + np.array(dlm4.result.filteredObs) - np.array(dlm5.result.filteredObs) + ), + 0.0, + ) def testPopout(self): dlm4 = dlm(self.data) - dlm4 + self.trend1 + dynamic(features = self.features, discount = 1) + dlm4 + self.trend1 + dynamic(features=self.features, discount=1) dlm4.fitForwardFilter() # the filtered step range should be (0, 19) self.assertEqual(dlm4.result.filteredSteps, [0, 19]) @@ -134,88 +154,95 @@ def testPopout(self): self.assertEqual(dlm4.result.filteredSteps, [0, -1]) dlm4.fitForwardFilter() - dlm5 = dlm(self.data[1 : 20]) - dlm5 + self.trend1 + dynamic(features = self.features[1 : 20], - discount = 1) + dlm5 = dlm(self.data[1:20]) + dlm5 + self.trend1 + dynamic(features=self.features[1:20], discount=1) dlm5.fitForwardFilter() # The two chain should have the same filtered obs - self.assertAlmostEqual(np.sum(np.array(dlm4.result.filteredObs) - - np.array(dlm5.result.filteredObs)), 0.0) - + self.assertAlmostEqual( + np.sum( + np.array(dlm4.result.filteredObs) - np.array(dlm5.result.filteredObs) + ), + 0.0, + ) def testAlter(self): dlm4 = dlm(self.data) - dlm4 + self.trend1 + dynamic(features = self.features, discount = 1) + dlm4 + self.trend1 + dynamic(features=self.features, discount=1) dlm4.fitForwardFilter() # the filtered step range should be (0, 19) self.assertEqual(dlm4.result.filteredSteps, [0, 19]) - dlm4.alter(date = 15, data = 1, component = 'main') + dlm4.alter(date=15, data=1, component="main") self.assertEqual(dlm4.result.filteredSteps, [0, 14]) dlm4.fitForwardFilter() newData = [0] * 9 + [1] + [0] * 10 newData[15] = 1 dlm5 = dlm(newData) - dlm5 + self.trend1 + dynamic(features = self.features, discount = 1) + dlm5 + self.trend1 + dynamic(features=self.features, discount=1) dlm5.fitForwardFilter() # The two chain should have the same filtered obs - self.assertAlmostEqual(np.sum(np.array(dlm4.result.filteredObs) - \ - np.array(dlm5.result.filteredObs)), 0.0) + self.assertAlmostEqual( + np.sum( + np.array(dlm4.result.filteredObs) - np.array(dlm5.result.filteredObs) + ), + 0.0, + ) # test alter the feature - dlm4.alter(date=0, data=[1,1], component='dynamic') - self.assertAlmostEqual(dlm4.builder.dynamicComponents['dynamic'].features[0], - [1, 1]) - + dlm4.alter(date=0, data=[1, 1], component="dynamic") + self.assertAlmostEqual( + dlm4.builder.dynamicComponents["dynamic"].features[0], [1, 1] + ) def testOneDayAheadPredictWithoutDynamic(self): self.dlm3.fitForwardFilter() (obs, var) = self.dlm3.predict(date=11) - self.assertAlmostEqual(obs, -6.0/7) - self.assertAlmostEqual(self.dlm3._predictModel.result.predictStatus, - [11, 12, [-6.0/7]]) + self.assertAlmostEqual(obs, -6.0 / 7) + self.assertAlmostEqual( + self.dlm3._predictModel.result.predictStatus, [11, 12, [-6.0 / 7]] + ) (obs, var) = self.dlm3.predict(date=2) - self.assertAlmostEqual(obs, 3.0/5) + self.assertAlmostEqual(obs, 3.0 / 5) # notice that the two latent states always sum up to 0 - self.assertAlmostEqual(self.dlm3._predictModel.result.predictStatus, - [2, 3, [3.0/5]]) - + self.assertAlmostEqual( + self.dlm3._predictModel.result.predictStatus, [2, 3, [3.0 / 5]] + ) def testOneDayAheadPredictWithDynamic(self): self.dlm4.fitForwardFilter() - featureDict = {'dynamic': 2.0} - (obs, var) = self.dlm4.predict(date=9, - featureDict=featureDict) - self.assertAlmostEqual(obs, 5.0/6 * 2) - + featureDict = {"dynamic": 2.0} + (obs, var) = self.dlm4.predict(date=9, featureDict=featureDict) + self.assertAlmostEqual(obs, 5.0 / 6 * 2) def testContinuePredictWithoutDynamic(self): self.dlm3.fitForwardFilter() (obs, var) = self.dlm3.predict(date=11) - self.assertAlmostEqual(self.dlm3._predictModel.result.predictStatus, - [11, 12, [-6.0/7]]) + self.assertAlmostEqual( + self.dlm3._predictModel.result.predictStatus, [11, 12, [-6.0 / 7]] + ) (obs, var) = self.dlm3.continuePredict() - self.assertAlmostEqual(self.dlm3._predictModel.result.predictStatus, - [11, 13, [-6.0/7, 6.0/7]]) - + self.assertAlmostEqual( + self.dlm3._predictModel.result.predictStatus, [11, 13, [-6.0 / 7, 6.0 / 7]] + ) def testContinuePredictWithDynamic(self): self.dlm4.fitForwardFilter() - featureDict = {'dynamic': [2.0]} - (obs, var) = self.dlm4.predict(date=9, - featureDict=featureDict) - self.assertAlmostEqual(self.dlm4._predictModel.result.predictStatus, - [9, 10, [5.0/6 * 2]]) + featureDict = {"dynamic": [2.0]} + (obs, var) = self.dlm4.predict(date=9, featureDict=featureDict) + self.assertAlmostEqual( + self.dlm4._predictModel.result.predictStatus, [9, 10, [5.0 / 6 * 2]] + ) - featureDict = {'dynamic': [3.0]} + featureDict = {"dynamic": [3.0]} (obs, var) = self.dlm4.continuePredict(featureDict=featureDict) - self.assertAlmostEqual(self.dlm4._predictModel.result.predictStatus, - [9, 11, [5.0/6 * 2, 5.0/6 * 3]]) - + self.assertAlmostEqual( + self.dlm4._predictModel.result.predictStatus, + [9, 11, [5.0 / 6 * 2, 5.0 / 6 * 3]], + ) def testPredictWithAutoReg(self): self.dlm5.fitForwardFilter() @@ -224,7 +251,6 @@ def testPredictWithAutoReg(self): (obs, var) = self.dlm5.continuePredict() self.assertAlmostEqual(obs[0, 0], 101.07480945) - def testPredictWithAutoReg2(self): self.dlm6.fitForwardFilter() (obs, var) = self.dlm6.predict(date=99) @@ -234,22 +260,21 @@ def testPredictWithAutoReg2(self): (obs, var) = self.dlm6.continuePredict() self.assertAlmostEqual(obs[0, 0], 102.0946503) - def testPredictNWithoutDynamic(self): self.dlm3.fitForwardFilter() (obs, var) = self.dlm3.predictN(N=2, date=11) - self.assertAlmostEqual(self.dlm3._predictModel.result.predictStatus, - [11, 13, [-6.0/7, 6.0/7]]) - + self.assertAlmostEqual( + self.dlm3._predictModel.result.predictStatus, [11, 13, [-6.0 / 7, 6.0 / 7]] + ) def testPredictNWithDynamic(self): self.dlm4.fitForwardFilter() - featureDict = {'dynamic': [[2.0], [3.0]]} - (obs, var) = self.dlm4.predictN(N=2, date=9, - featureDict=featureDict) - self.assertAlmostEqual(self.dlm4._predictModel.result.predictStatus, - [9, 11, [5.0/6 * 2, 5.0/6 * 3]]) - + featureDict = {"dynamic": [[2.0], [3.0]]} + (obs, var) = self.dlm4.predictN(N=2, date=9, featureDict=featureDict) + self.assertAlmostEqual( + self.dlm4._predictModel.result.predictStatus, + [9, 11, [5.0 / 6 * 2, 5.0 / 6 * 3]], + ) def testPredictNWithAutoReg(self): self.dlm5.fitForwardFilter() @@ -257,93 +282,85 @@ def testPredictNWithAutoReg(self): self.assertAlmostEqual(obs[0], 100.03682874) self.assertAlmostEqual(obs[1], 101.07480945) - def testPredictNWithDynamicMatrixInput(self): self.dlm4.fitForwardFilter() - featureDict = {'dynamic': np.array([[2.0], [3.0]])} - (obs, var) = self.dlm4.predictN(N=2, date=9, - featureDict=featureDict) - self.assertAlmostEqual(self.dlm4._predictModel.result.predictStatus, - [9, 11, [5.0/6 * 2, 5.0/6 * 3]]) - + featureDict = {"dynamic": np.array([[2.0], [3.0]])} + (obs, var) = self.dlm4.predictN(N=2, date=9, featureDict=featureDict) + self.assertAlmostEqual( + self.dlm4._predictModel.result.predictStatus, + [9, 11, [5.0 / 6 * 2, 5.0 / 6 * 3]], + ) def testPredictionNotChangeModel(self): timeSeries = [1, 2, 1, 5, 3, 5, 4, 8, 1, 2] dlm1 = dlm(timeSeries) + trend(degree=2, discount=0.95) dlm1.fitForwardFilter() - (obs1, var1) = dlm1.predictN(N=1, date=dlm1.n-1) + (obs1, var1) = dlm1.predictN(N=1, date=dlm1.n - 1) dlm2 = dlm([]) + trend(degree=2, discount=0.95) for d in timeSeries: - dlm2.append([d], component='main') + dlm2.append([d], component="main") dlm2.fitForwardFilter() - (obs2, var2) = dlm2.predictN(N=1, date=dlm2.n-1) - + (obs2, var2) = dlm2.predictN(N=1, date=dlm2.n - 1) + self.assertAlmostEqual(obs1, obs2) self.assertAlmostEqual(var1, var2) - def testGetLatentState(self): # for forward filter self.dlm5.fitForwardFilter() filteredTrend = self.dlm5.getLatentState( - filterType='forwardFilter', name='trend') + filterType="forwardFilter", name="trend" + ) diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(filteredTrend[i][0] - - self.dlm5.result.filteredState[i][0, 0]) + diff += abs(filteredTrend[i][0] - self.dlm5.result.filteredState[i][0, 0]) self.assertAlmostEqual(diff, 0) # for backward smoother self.dlm5.fitBackwardSmoother() smoothedTrend = self.dlm5.getLatentState( - filterType='backwardSmoother', name='trend') + filterType="backwardSmoother", name="trend" + ) diff = 0.0 for i in range(len(smoothedTrend)): - diff += abs(smoothedTrend[i][0] - - self.dlm5.result.smoothedState[i][0, 0]) + diff += abs(smoothedTrend[i][0] - self.dlm5.result.smoothedState[i][0, 0]) self.assertAlmostEqual(diff, 0) - def testGetLatentCov(self): # for forward filter self.dlm5.fitForwardFilter() - filteredTrend = self.dlm5.getLatentCov( - filterType='forwardFilter', name='trend') + filteredTrend = self.dlm5.getLatentCov(filterType="forwardFilter", name="trend") diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(filteredTrend[i][0, 0] - - self.dlm5.result.filteredCov[i][0, 0]) + diff += abs(filteredTrend[i][0, 0] - self.dlm5.result.filteredCov[i][0, 0]) self.assertAlmostEqual(diff, 0) # for backward smoother self.dlm5.fitBackwardSmoother() smoothedTrend = self.dlm5.getLatentCov( - filterType='backwardSmoother', name='trend') + filterType="backwardSmoother", name="trend" + ) diff = 0.0 for i in range(len(smoothedTrend)): - diff += abs(smoothedTrend[i][0, 0] - - self.dlm5.result.smoothedCov[i][0, 0]) + diff += abs(smoothedTrend[i][0, 0] - self.dlm5.result.smoothedCov[i][0, 0]) self.assertAlmostEqual(diff, 0) - def testGetMean(self): # for forward filter self.dlm5.fitForwardFilter() - filteredTrend = self.dlm5.getMean(filterType='forwardFilter') + filteredTrend = self.dlm5.getMean(filterType="forwardFilter") self.assertEqual(len(filteredTrend), self.dlm5.n) diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(filteredTrend[i] - - self.dlm5.result.filteredObs[i][0, 0]) + diff += abs(filteredTrend[i] - self.dlm5.result.filteredObs[i][0, 0]) self.assertAlmostEqual(diff, 0) # for component with forward filter - arTrend = self.dlm5.getMean(filterType='forwardFilter', - name='ar2') + arTrend = self.dlm5.getMean(filterType="forwardFilter", name="ar2") trueAr = [item[1, 0] for item in self.dlm5.result.filteredState] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] @@ -355,18 +372,16 @@ def testGetMean(self): # for backward smoother self.dlm5.fitBackwardSmoother() - filteredTrend = self.dlm5.getMean(filterType='backwardSmoother') + filteredTrend = self.dlm5.getMean(filterType="backwardSmoother") diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(filteredTrend[i] - - self.dlm5.result.smoothedObs[i][0, 0]) + diff += abs(filteredTrend[i] - self.dlm5.result.smoothedObs[i][0, 0]) self.assertAlmostEqual(diff, 0) # for component with backward smoother - arTrend = self.dlm5.getMean(filterType='backwardSmoother', - name='ar2') + arTrend = self.dlm5.getMean(filterType="backwardSmoother", name="ar2") trueAr = [item[1, 0] for item in self.dlm5.result.smoothedState] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] @@ -376,23 +391,20 @@ def testGetMean(self): diff += abs(arTrend[i] - trueAr[i]) self.assertAlmostEqual(diff, 0) - def testGetVar(self): # for forward filter self.dlm5.fitForwardFilter() - filteredTrend = self.dlm5.getVar(filterType='forwardFilter') + filteredTrend = self.dlm5.getVar(filterType="forwardFilter") self.assertEqual(len(filteredTrend), self.dlm5.n) diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(filteredTrend[i] - - self.dlm5.result.filteredObsVar[i][0, 0]) + diff += abs(filteredTrend[i] - self.dlm5.result.filteredObsVar[i][0, 0]) self.assertAlmostEqual(diff, 0) # for component with forward filter - arTrend = self.dlm5.getVar(filterType='forwardFilter', - name='ar2') + arTrend = self.dlm5.getVar(filterType="forwardFilter", name="ar2") trueAr = [item[1, 1] for item in self.dlm5.result.filteredCov] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] * comp.evaluation.T @@ -404,18 +416,16 @@ def testGetVar(self): # for backward smoother self.dlm5.fitBackwardSmoother() - filteredTrend = self.dlm5.getVar(filterType='backwardSmoother') + filteredTrend = self.dlm5.getVar(filterType="backwardSmoother") diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(filteredTrend[i] - - self.dlm5.result.smoothedObsVar[i][0, 0]) + diff += abs(filteredTrend[i] - self.dlm5.result.smoothedObsVar[i][0, 0]) self.assertAlmostEqual(diff, 0) # for component with backward smoother - arTrend = self.dlm5.getVar(filterType='backwardSmoother', - name='ar2') + arTrend = self.dlm5.getVar(filterType="backwardSmoother", name="ar2") trueAr = [item[1, 1] for item in self.dlm5.result.smoothedCov] - comp = self.dlm5.builder.automaticComponents['ar2'] + comp = self.dlm5.builder.automaticComponents["ar2"] for i in range(len(trueAr)): comp.updateEvaluation(i, self.data5) trueAr[i] = comp.evaluation * trueAr[i] * comp.evaluation.T @@ -425,15 +435,13 @@ def testGetVar(self): diff += abs(arTrend[i] - trueAr[i]) self.assertAlmostEqual(diff, 0) - def testGetMSE(self): self.dlm1.stableMode(False) self.dlm1.fitForwardFilter() mse1 = self.dlm1.getMSE() mse_expect = 0 for i in range(20): - mse_expect += (self.dlm1.result.predictedObs[i] - - self.data[i]) ** 2 + mse_expect += (self.dlm1.result.predictedObs[i] - self.data[i]) ** 2 mse_expect /= 20 self.assertAlmostEqual(mse1, mse_expect) @@ -441,35 +449,31 @@ def testGetMSE(self): self.dlm2.fitForwardFilter() self.dlm2.result.filteredSteps = (0, 19) mse2 = self.dlm2.getMSE() - mse_expect = 2.0/20 + mse_expect = 2.0 / 20 self.assertAlmostEqual(mse2, mse_expect) - def testGetResidual(self): # for forward filter - filter_type = 'forwardFilter' + filter_type = "forwardFilter" self.dlm5.fitForwardFilter() filteredTrend = self.dlm5.getMean(filterType=filter_type) filteredResidual = self.dlm5.getResidual(filterType=filter_type) self.assertEqual(len(filteredResidual), self.dlm5.n) diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(- filteredTrend[i] - filteredResidual[i] + - self.dlm5.data[i]) + diff += abs(-filteredTrend[i] - filteredResidual[i] + self.dlm5.data[i]) self.assertAlmostEqual(diff, 0) # for backward smoother - filter_type = 'backwardSmoother' + filter_type = "backwardSmoother" self.dlm5.fitBackwardSmoother() filteredTrend = self.dlm5.getMean(filterType=filter_type) filteredResidual = self.dlm5.getResidual(filterType=filter_type) diff = 0.0 for i in range(len(filteredTrend)): - diff += abs(- filteredTrend[i] - filteredResidual[i] + - self.dlm5.data[i]) - self.assertAlmostEqual(diff, 0) - + diff += abs(-filteredTrend[i] - filteredResidual[i] + self.dlm5.data[i]) + self.assertAlmostEqual(diff, 0) def testTune(self): # just make sure the tune can run @@ -479,5 +483,6 @@ def testTune(self): def testBuildFromBuilder(self): self.dlm6.fit() -if __name__ == '__main__': + +if __name__ == "__main__": unittest.main() diff --git a/tests/tuner/testDlmTuner.py b/tests/tuner/testDlmTuner.py index 40e562a7..3ce6c515 100644 --- a/tests/tuner/testDlmTuner.py +++ b/tests/tuner/testDlmTuner.py @@ -8,27 +8,26 @@ class testModelTuner(unittest.TestCase): - - def setUp(self): self.mydlm = dlm(np.random.random(100)) + trend(2, discount=0.95) self.mytuner = modelTuner() - def testFind_gradient(self): mydlm2 = deepcopy(self.mydlm) self.mydlm.fitForwardFilter() mydlm2.fitForwardFilter() mse0 = mydlm2._getMSE() - mydlm2._setDiscounts(list(map(lambda x: x + self.mytuner.err, - self.mydlm._getDiscounts()))) + mydlm2._setDiscounts( + list(map(lambda x: x + self.mytuner.err, self.mydlm._getDiscounts())) + ) mydlm2.fitForwardFilter() mse1 = mydlm2._getMSE() expect_gradient = (mse1 - mse0) / self.mytuner.err self.assertAlmostEqual( - expect_gradient, self.mytuner.find_gradient( - self.mydlm._getDiscounts(), self.mydlm)) + expect_gradient, + self.mytuner.find_gradient(self.mydlm._getDiscounts(), self.mydlm), + ) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main() diff --git a/tests/tuner/test_dlmTune.py b/tests/tuner/test_dlmTune.py index 8627cc54..703b0f2f 100644 --- a/tests/tuner/test_dlmTune.py +++ b/tests/tuner/test_dlmTune.py @@ -7,8 +7,6 @@ class test_dlmTune(unittest.TestCase): - - def setUp(self): self.data = [0] * 9 + [1] + [0] * 10 self.data5 = range(100) @@ -18,57 +16,57 @@ def setUp(self): self.dlm7 = _dlmTune([0, 1, None, 1, 0, 1, -1]) self.dlm1.builder + trend(degree=0, discount=1, w=1.0) self.dlm2.builder + trend(degree=0, discount=1e-12, w=1.0) - self.dlm6.builder + trend(degree=0, discount=0.9, w=1.0) + \ - seasonality(period=2, discount=0.8, w=1.0) + \ - autoReg(degree=3, discount=1.0) + ( + self.dlm6.builder + + trend(degree=0, discount=0.9, w=1.0) + + seasonality(period=2, discount=0.8, w=1.0) + + autoReg(degree=3, discount=1.0) + ) self.dlm7.builder + trend(degree=0, discount=1, w=1.0) self.dlm1._initialize() self.dlm2._initialize() self.dlm6._initialize() self.dlm7._initialize() - self.dlm1.options.innovationType='whole' - self.dlm2.options.innovationType='whole' - self.dlm6.options.innovationType='whole' - self.dlm7.options.innovationType='whole' - + self.dlm1.options.innovationType = "whole" + self.dlm2.options.innovationType = "whole" + self.dlm6.options.innovationType = "whole" + self.dlm7.options.innovationType = "whole" def testComputeMSE(self): self.dlm1._forwardFilter(start=0, end=19, renew=False) - self.dlm1.result.filteredSteps=(0, 19) + self.dlm1.result.filteredSteps = (0, 19) mse1 = self.dlm1._getMSE() mse_expect = 0 for i in range(20): - mse_expect += (self.dlm1.result.predictedObs[i] - - self.data[i]) ** 2 + mse_expect += (self.dlm1.result.predictedObs[i] - self.data[i]) ** 2 mse_expect /= 20 self.assertAlmostEqual(mse1, mse_expect) self.dlm2._forwardFilter(start=0, end=19, renew=False) - self.dlm2.result.filteredSteps=(0, 19) + self.dlm2.result.filteredSteps = (0, 19) mse2 = self.dlm2._getMSE() - mse_expect = 2.0/20 + mse_expect = 2.0 / 20 self.assertAlmostEqual(mse2, mse_expect) # Test missing data self.dlm7._forwardFilter(start=0, end=6, renew=False) - self.dlm7.result.filteredSteps=(0, 6) + self.dlm7.result.filteredSteps = (0, 6) mse3 = self.dlm7._getMSE() mse_expect = 0 for i in range(7): if self.dlm7.data[i] is not None: - mse_expect += (self.dlm7.result.predictedObs[i] - - self.dlm7.data[i]) ** 2 + mse_expect += ( + self.dlm7.result.predictedObs[i] - self.dlm7.data[i] + ) ** 2 mse_expect /= 7 self.assertAlmostEqual(mse3, mse_expect) - def testGetDiscount(self): discounts = self.dlm6._getDiscounts() self.assertTrue(0.9 in discounts) self.assertTrue(0.8 in discounts) self.assertTrue(1.0 in discounts) - def testSetDiscount(self): self.dlm6._setDiscounts([0.0, 0.1, 0.2], False) self.assertTrue(0.0 in self.dlm6.builder.discount) @@ -78,13 +76,15 @@ def testSetDiscount(self): self.assertTrue(0.8 not in self.dlm6.builder.discount) self.assertTrue(1.0 not in self.dlm6.builder.discount) - self.assertAlmostEqual(self.dlm6.builder.staticComponents['trend'].discount, - 0.9) + self.assertAlmostEqual( + self.dlm6.builder.staticComponents["trend"].discount, 0.9 + ) self.dlm6._setDiscounts([0.0, 0.1, 0.2], True) - self.assertTrue(self.dlm6.builder.staticComponents['trend'].discount in - [0.0, 0.1, 0.2]) + self.assertTrue( + self.dlm6.builder.staticComponents["trend"].discount in [0.0, 0.1, 0.2] + ) -if __name__ == '__main__': +if __name__ == "__main__": unittest.main()