Commit 1ddb618f authored by Andrew Quinn's avatar Andrew Quinn
Browse files

Merge branch 'mh1/estimators' into 'master'

Allow passing sklearn estimator to fit_model

See merge request analysis/sails!63
parents 5c129835 4155c617
Loading
Loading
Loading
Loading
Loading
+1 −1
Original line number Diff line number Diff line
@@ -9,7 +9,7 @@ build_test:
  script:
    - apt update
    - apt upgrade -y
    - apt install -y python3-setuptools python3-numpy python3-scipy python3-matplotlib python3-h5py python3-sphinx python3-sphinx-rtd-theme ipython3 flake8 python3-setuptools python3-pytest python3-pytest-cov python3-wheel codespell
    - apt install -y python3-setuptools python3-numpy python3-scipy python3-matplotlib python3-h5py python3-sklearn python3-sphinx python3-sphinx-rtd-theme ipython3 flake8 python3-setuptools python3-pytest python3-pytest-cov python3-wheel codespell
    - flake8
    - make spell
    - python3 setup.py build
+33 −13
Original line number Diff line number Diff line
@@ -25,15 +25,35 @@ We then fit a model using Ordinary Least Squared as implemented in SAILS.
   # Fit an autoregressive model with order 3
   sails_model = sails.OLSLinearModel.fit_model(X, np.arange(4))

Now we will create a new model fit class based on
``sails.modelfit.AbstractLinearModel``. This is a base class which contains a
number of methods and properties to store and compute information on a model.
The ``AbstractLinearModel`` is not usable on its own as the fit_model method is
not implemented. When classes such as ``OLSLinearModel`` are defined in SAILS,
they inherit from ``AbstractLinearModel`` to define the helper functions before a
specific ``fit_model`` method is defined. We can do the same to define a custom
model fit class using an external package. We will first create a new class
which inherits from ``AbstractLinearModel`` and then define a ``fit_model`` method
Before we continue, we should note that you can pass an existing ``sklearn``
estimator (for example, ``sklearn.linear_model.LinearRegression`` as the
``estimator`` parameter to the ``fit_model`` function of the ``OLSLinearModel``
class.  If you do this, you should not fit the intercept in the model.  For
instance:

.. code-block:: python

   import sklearn.linear_model

   # Fit an autoregressive model using SKLearn's LinearRegressor
   estimator = sklearn.linear_model.LinearRegression(fit_intercept=False)
   sails_model = sails.OLSLinearModel.fit_model(X, np.arange(4), estimator=estimator)

The above will give the same answers as the default method (which calculates
the parameters using the normal equations).  You can, however, extend this
approach to use, for instance, ridge or lasso-regression using the relevant
classes (`sklearn.linear_model.Ridge` or `sklearn.linear_model.Lasso`).

To go beyond what is available using the previous options, we can create a new
model fit class based on ``sails.modelfit.AbstractLinearModel``. This is a base
class which contains a number of methods and properties to store and compute
information on a model.  The ``AbstractLinearModel`` is not usable on its own
as the fit_model method is not implemented. When classes such as
``OLSLinearModel`` are defined in SAILS, they inherit from
``AbstractLinearModel`` to define the helper functions before a specific
``fit_model`` method is defined. We can do the same to define a custom model
fit class using an external package. We will first create a new class which
inherits from ``AbstractLinearModel`` and then define a ``fit_model`` method
which computes the model fit and stores the outputs in a standard form.

Here is our custom model fit class, each section is described in the comments in the code.
+1 −0
Original line number Diff line number Diff line
@@ -2,4 +2,5 @@ numpy>=1.16
scipy>=1.1.0
matplotlib>=3.0.2
h5py>=2.8.0
sklearn>=0.20
sphinx_rtd_theme
+13 −3
Original line number Diff line number Diff line
@@ -423,7 +423,7 @@ class OLSLinearModel(AbstractLinearModel):
    """A class implementing ordinary least squares linear model fit"""

    @classmethod
    def fit_model(cls, data, delay_vect):
    def fit_model(cls, data, delay_vect, estimator=None):
        """This is a class method which fits a linear model and returns a
        populated :class:`OLSLinearModel` instance containing the fitted model.

@@ -433,6 +433,9 @@ class OLSLinearModel(AbstractLinearModel):
            The data to be used to compute the model fit
        delay_vect : ndarray
            A vector of lag indices defining the lags to fit
        estimator : None or sklearn class
            If None, fit using standard OLS normal equations.
            If set to an appropriate sklearn class, use that estimator.


        Returns
@@ -467,8 +470,15 @@ class OLSLinearModel(AbstractLinearModel):
        # [ A[1, 1, 1] A[1, 1, 2] ... A[1, 1, p] A[1, 2, 1] A[1, 2, 2] ... A[1, 2, p] ]
        # [ A[2, 1, 1] A[2, 1, 2] ... A[2, 1, p] A[2, 2, 1] A[2, 2, 2] ... A[2, 2, p] ]

        # Normal equations
        if estimator is None:
            # Use normal equations
            B = np.linalg.inv(X.T.dot(X)).dot(X.T).dot(Y)
        else:
            # Fit the passed in estimator with our data
            estimator.fit(X, Y)

            # The coefficients are transposed in the sklearn implementation
            B = estimator.coef_.T

        # Reshape output
        parameters = np.zeros((nsignals, nsignals, maxorder+1))
+20 −2
Original line number Diff line number Diff line
@@ -13,15 +13,18 @@ __all__ = []

class test_model_fits(unittest.TestCase):

    def run_fit_checks(self, freq, fitfunc):
    def run_fit_checks(self, freq, fitfunc, extraargs=None):
        from ..modal import MvarModalDecomposition

        if extraargs is None:
            extraargs = {}

        sample_rate = 100
        seconds = 50
        time_vect = np.linspace(0, seconds, sample_rate * seconds)
        x = np.sin(2*np.pi*freq*time_vect)

        model = fitfunc.fit_model(x[None, :, None], np.arange(3))
        model = fitfunc.fit_model(x[None, :, None], np.arange(3), **extraargs)
        modes = MvarModalDecomposition.initialise(model, sample_rate=sample_rate)

        # Compute analytic modes
@@ -49,6 +52,21 @@ class test_model_fits(unittest.TestCase):
        self.run_fit_checks(10, fit)
        self.run_fit_checks(40, fit)

    def test_ols_sklearn(self):
        """Test the Ordinary Least Squares model fits using sklearn"""

        import sklearn.linear_model

        from ..modelfit import OLSLinearModel

        fit = OLSLinearModel
        estimator = sklearn.linear_model.LinearRegression(fit_intercept=False)
        extraargs = {'estimator': estimator}

        self.run_fit_checks(24, fit, extraargs=extraargs)
        self.run_fit_checks(10, fit, extraargs=extraargs)
        self.run_fit_checks(40, fit, extraargs=extraargs)

    def test_vieiramorf(self):
        """Test the Vieira-Morf model fits"""