add subclass for levenberg-marquardt method. update initial search sp…

…ace. bump version. (#102)

.gitignore

@@ -5,6 +5,8 @@ __pycache__/
 .idea
 .DS_Store
 
+*/.ipynb_checkpoints/
+
 example.py
 
 outcmaes

lppls/lppls.py

@@ -39,7 +39,6 @@ def func_restricted(self, x, *args):
         Returns:
             (float)
         """
-
         tc = x[0]
         m = x[1]
         w = x[2]
@@ -53,7 +52,6 @@ def func_restricted(self, x, *args):
         delta = [self.lppls(t, tc, m, w, a, b, c1, c2) for t in observations[0, :]]
         delta = np.subtract(delta, observations[1, :])
         delta = np.power(delta, 2)
-
         return np.sum(delta)
 
     @staticmethod
@@ -116,22 +114,20 @@ def fit(self, max_searches, minimizer='Nelder-Mead', obs=None):
         if obs is None:
             obs = self.observations
 
-
-        # print('obs',obs)
         search_count = 0
         # find bubble
         while search_count < max_searches:
-            tc_init_min, tc_init_max = self._get_tc_bounds(obs, 0.50, 0.50)
+            # tc_init_min, tc_init_max = self._get_tc_bounds(obs, 0.50, 0.50)
             t1 = obs[0, 0]
             t2 = obs[0, -1]
 
             # @TODO make configurable
             # set random initialization limits for non-linear params
             init_limits = [
-                (max(t2 - 60, t2 - 0.5 * (t2 - t1)), min(t2 + 252, t2 + 0.5 * (t2 - t1))),  # tc
                 # (tc_init_min, tc_init_max),
-                (0.0, 1.0),  # m
-                (2.0, 15.0),  # ω
+                (t2 - 0.2 * (t2 - t1), t2 + 0.2 * (t2 - t1)),  # tc
+                (0.1, 1.0),  # m
+                (6.0, 13.0),  # ω
             ]
 
             # randomly choose vals within bounds for non-linear params

lppls/lppls_lm.py

@@ -0,0 +1,48 @@
+import numpy as np
+from scipy.optimize import least_squares
+from lppls.lppls import LPPLS
+
+
+class LPPLS_LM(LPPLS):
+
+    def func_restricted(self, x, obs):
+        tc = x[0]
+        m = x[1]
+        w = x[2]
+
+        rM = self.matrix_equation(obs, tc, m, w)
+        a, b, c1, c2 = rM[:, 0].tolist()
+
+        delta = [self.lppls(t, tc, m, w, a, b, c1, c2) for t in obs[0, :]]
+        residuals = np.subtract(delta, obs[1, :])
+        return residuals  # return the array of residuals
+
+    def estimate_params(self, observations, seed, minimizer=None):
+        """
+        Overrides the estimate_params method to use least_squares with 'lm' method.
+        Args:
+            observations (np.ndarray): The observed time-series data.
+            seed (list): Initial guess for time-critical, omega, and m.
+        Returns:
+            tc, m, w, a, b, c, c1, c2
+        """
+        # Define a wrapper function for least_squares
+        def wrapper(x):
+            return self.func_restricted(x, observations)
+
+        # Use least_squares with the Levenberg-Marquardt method
+        result = least_squares(wrapper, seed, method='lm')
+
+        if result.success:
+            tc, m, w = result.x
+            rM = self.matrix_equation(observations, tc, m, w)
+            a, b, c1, c2 = rM[:, 0].tolist()
+            c = self.get_c(c1, c2)
+
+            # Store fitted parameters
+            for coef in ['tc', 'm', 'w', 'a', 'b', 'c', 'c1', 'c2']:
+                self.coef_[coef] = eval(coef)
+
+            return tc, m, w, a, b, c, c1, c2
+        else:
+            raise ValueError("Parameter estimation failed.")

setup.py

@@ -4,7 +4,7 @@
     long_description = fh.read()
 
 setuptools.setup(name='lppls',
-      version='0.6.10',
+      version='0.6.11',
       description='A Python module for fitting the LPPLS model to data.',
       packages=['lppls'],
       author='Josh Nielsen',