|
2 | 2 | Fitting in EasyScience |
3 | 3 | ====================== |
4 | 4 |
|
| 5 | +EasyScience provides a flexible and powerful fitting framework that supports multiple optimization backends. |
| 6 | +This guide covers both basic usage for users wanting to fit their data, and advanced patterns for developers building scientific components. |
| 7 | + |
| 8 | +Overview |
| 9 | +-------- |
| 10 | + |
| 11 | +The EasyScience fitting system consists of: |
| 12 | + |
| 13 | +* **Parameters**: Scientific values with units, bounds, and fitting capabilities |
| 14 | +* **Models**: Objects containing parameters, inheriting from ``ObjBase`` |
| 15 | +* **Fitter**: The main fitting engine supporting multiple minimizers |
| 16 | +* **Minimizers**: Backend optimization engines (LMFit, Bumps, DFO-LS) |
| 17 | + |
| 18 | +Quick Start |
| 19 | +----------- |
| 20 | + |
| 21 | +Basic Parameter and Model Setup |
| 22 | +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
| 23 | + |
| 24 | +.. code-block:: python |
| 25 | +
|
| 26 | + import numpy as np |
| 27 | + from easyscience import ObjBase, Parameter, Fitter |
| 28 | +
|
| 29 | + # Create a simple model with fittable parameters |
| 30 | + class SineModel(ObjBase): |
| 31 | + def __init__(self, amplitude_val=1.0, frequency_val=1.0, phase_val=0.0): |
| 32 | + amplitude = Parameter("amplitude", amplitude_val, min=0, max=10) |
| 33 | + frequency = Parameter("frequency", frequency_val, min=0.1, max=5) |
| 34 | + phase = Parameter("phase", phase_val, min=-np.pi, max=np.pi) |
| 35 | + super().__init__("sine_model", amplitude=amplitude, frequency=frequency, phase=phase) |
| 36 | + |
| 37 | + def __call__(self, x): |
| 38 | + return self.amplitude.value * np.sin(2 * np.pi * self.frequency.value * x + self.phase.value) |
| 39 | +
|
| 40 | +Basic Fitting Example |
| 41 | +~~~~~~~~~~~~~~~~~~~~~ |
| 42 | + |
| 43 | +.. code-block:: python |
| 44 | +
|
| 45 | + # Create test data |
| 46 | + x_data = np.linspace(0, 2, 100) |
| 47 | + true_model = SineModel(amplitude_val=2.5, frequency_val=1.5, phase_val=0.5) |
| 48 | + y_data = true_model(x_data) + 0.1 * np.random.normal(size=len(x_data)) |
| 49 | + |
| 50 | + # Create model to fit with initial guesses |
| 51 | + fit_model = SineModel(amplitude_val=1.0, frequency_val=1.0, phase_val=0.0) |
| 52 | + |
| 53 | + # Set which parameters to fit (unfix them) |
| 54 | + fit_model.amplitude.fixed = False |
| 55 | + fit_model.frequency.fixed = False |
| 56 | + fit_model.phase.fixed = False |
| 57 | + |
| 58 | + # Create fitter and perform fit |
| 59 | + fitter = Fitter(fit_model, fit_model) |
| 60 | + result = fitter.fit(x=x_data, y=y_data) |
| 61 | + |
| 62 | + # Access results |
| 63 | + print(f"Chi-squared: {result.chi2}") |
| 64 | + print(f"Fitted amplitude: {fit_model.amplitude.value} ± {fit_model.amplitude.error}") |
| 65 | + print(f"Fitted frequency: {fit_model.frequency.value} ± {fit_model.frequency.error}") |
| 66 | +
|
| 67 | +Available Minimizers |
| 68 | +------------------- |
| 69 | + |
| 70 | +EasyScience supports multiple optimization backends: |
| 71 | + |
| 72 | +.. code-block:: python |
| 73 | +
|
| 74 | + from easyscience import AvailableMinimizers |
| 75 | + |
| 76 | + # View all available minimizers |
| 77 | + fitter = Fitter(model, model) |
| 78 | + print(fitter.available_minimizers) |
| 79 | + # Output: ['LMFit', 'LMFit_leastsq', 'LMFit_powell', 'Bumps', 'Bumps_simplex', 'DFO', 'DFO_leastsq'] |
| 80 | +
|
| 81 | +Switching Minimizers |
| 82 | +~~~~~~~~~~~~~~~~~~~ |
| 83 | + |
| 84 | +.. code-block:: python |
| 85 | +
|
| 86 | + # Use LMFit (default) |
| 87 | + fitter.switch_minimizer(AvailableMinimizers.LMFit) |
| 88 | + result1 = fitter.fit(x=x_data, y=y_data) |
| 89 | + |
| 90 | + # Switch to Bumps |
| 91 | + fitter.switch_minimizer(AvailableMinimizers.Bumps) |
| 92 | + result2 = fitter.fit(x=x_data, y=y_data) |
| 93 | + |
| 94 | + # Use DFO for derivative-free optimization |
| 95 | + fitter.switch_minimizer(AvailableMinimizers.DFO) |
| 96 | + result3 = fitter.fit(x=x_data, y=y_data) |
| 97 | +
|
| 98 | +Parameter Management |
| 99 | +------------------- |
| 100 | + |
| 101 | +Setting Bounds and Constraints |
| 102 | +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
| 103 | + |
| 104 | +.. code-block:: python |
| 105 | +
|
| 106 | + # Parameter with bounds |
| 107 | + param = Parameter(name="amplitude", value=1.0, min=0.0, max=10.0, unit="m") |
| 108 | + |
| 109 | + # Fix parameter (exclude from fitting) |
| 110 | + param.fixed = True |
| 111 | + |
| 112 | + # Unfix parameter (include in fitting) |
| 113 | + param.fixed = False |
| 114 | + |
| 115 | + # Change bounds dynamically |
| 116 | + param.min = 0.5 |
| 117 | + param.max = 8.0 |
| 118 | +
|
| 119 | +Parameter Dependencies |
| 120 | +~~~~~~~~~~~~~~~~~~~~~ |
| 121 | + |
| 122 | +Parameters can depend on other parameters through expressions: |
| 123 | + |
| 124 | +.. code-block:: python |
| 125 | +
|
| 126 | + # Create independent parameters |
| 127 | + length = Parameter("length", 10.0, unit="m", min=1, max=100) |
| 128 | + width = Parameter("width", 5.0, unit="m", min=1, max=50) |
| 129 | + |
| 130 | + # Create dependent parameter |
| 131 | + area = Parameter.from_dependency( |
| 132 | + name="area", |
| 133 | + dependency_expression="length * width", |
| 134 | + dependency_map={"length": length, "width": width} |
| 135 | + ) |
| 136 | + |
| 137 | + # When length or width changes, area updates automatically |
| 138 | + length.value = 15.0 |
| 139 | + print(area.value) # Will be 75.0 (15 * 5) |
| 140 | +
|
| 141 | +Using make_dependent_on() Method |
| 142 | +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
| 143 | + |
| 144 | +You can also make an existing parameter dependent on other parameters using the ``make_dependent_on()`` method. This is useful when you want to convert an independent parameter into a dependent one: |
| 145 | + |
| 146 | +.. code-block:: python |
| 147 | +
|
| 148 | + # Create independent parameters |
| 149 | + radius = Parameter("radius", 5.0, unit="m", min=1, max=20) |
| 150 | + height = Parameter("height", 10.0, unit="m", min=1, max=50) |
| 151 | + volume = Parameter("volume", 100.0, unit="m³") # Initially independent |
| 152 | + pi = Parameter("pi", 3.14159, fixed=True) # Constant parameter |
| 153 | +
|
| 154 | + # Make volume dependent on radius and height |
| 155 | + volume.make_dependent_on( |
| 156 | + dependency_expression="pi * radius**2 * height", |
| 157 | + dependency_map={"radius": radius, "height": height, "pi": pi} |
| 158 | + ) |
| 159 | +
|
| 160 | + # Now volume automatically updates when radius or height changes |
| 161 | + radius.value = 8.0 |
| 162 | + print(f"New volume: {volume.value:.2f} m³") # Automatically calculated |
| 163 | +
|
| 164 | + # The parameter becomes dependent and cannot be set directly |
| 165 | + try: |
| 166 | + volume.value = 200.0 # This will raise an AttributeError |
| 167 | + except AttributeError: |
| 168 | + print("Cannot set value of dependent parameter directly") |
| 169 | +
|
| 170 | +**What to expect:** |
| 171 | + |
| 172 | +- The parameter becomes **dependent** and its ``independent`` property becomes ``False`` |
| 173 | +- You **cannot directly set** the value, bounds, or variance of a dependent parameter |
| 174 | +- The parameter's value is **automatically recalculated** whenever any of its dependencies change |
| 175 | +- Dependent parameters **cannot be fitted** (they are automatically fixed) |
| 176 | +- The original value, unit, variance, min, and max are **overwritten** by the dependency calculation |
| 177 | +- You can **revert to independence** using the ``make_independent()`` method if needed |
| 178 | + |
| 179 | +Advanced Fitting Options |
| 180 | +----------------------- |
| 181 | + |
| 182 | +Setting Tolerances and Limits |
| 183 | +~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ |
| 184 | + |
| 185 | +.. code-block:: python |
| 186 | +
|
| 187 | + fitter = Fitter(model, model) |
| 188 | + |
| 189 | + # Set convergence tolerance |
| 190 | + fitter.tolerance = 1e-8 |
| 191 | + |
| 192 | + # Limit maximum function evaluations |
| 193 | + fitter.max_evaluations = 1000 |
| 194 | + |
| 195 | + # Perform fit with custom settings |
| 196 | + result = fitter.fit(x=x_data, y=y_data) |
| 197 | +
|
| 198 | +Using Weights |
| 199 | +~~~~~~~~~~~~ |
| 200 | + |
| 201 | +.. code-block:: python |
| 202 | +
|
| 203 | + # Define weights (inverse variance) |
| 204 | + weights = 1.0 / errors**2 # where errors are your data uncertainties |
| 205 | + |
| 206 | + # Fit with weights |
| 207 | + result = fitter.fit(x=x_data, y=y_data, weights=weights) |
| 208 | +
|
| 209 | +Multidimensional Fitting |
| 210 | +~~~~~~~~~~~~~~~~~~~~~~~~ |
| 211 | + |
| 212 | +.. code-block:: python |
| 213 | +
|
| 214 | + class AbsSin2D(ObjBase): |
| 215 | + def __init__(self, offset_val=0.0, phase_val=0.0): |
| 216 | + offset = Parameter("offset", offset_val) |
| 217 | + phase = Parameter("phase", phase_val) |
| 218 | + super().__init__("sin2D", offset=offset, phase=phase) |
| 219 | + |
| 220 | + def __call__(self, x): |
| 221 | + X, Y = x[:, 0], x[:, 1] # x is 2D array |
| 222 | + return np.abs(np.sin(self.phase.value * X + self.offset.value)) * \ |
| 223 | + np.abs(np.sin(self.phase.value * Y + self.offset.value)) |
| 224 | + |
| 225 | + # Create 2D data |
| 226 | + x_2d = np.column_stack([x_grid.ravel(), y_grid.ravel()]) |
| 227 | + |
| 228 | + # Fit 2D model |
| 229 | + model_2d = AbsSin2D(offset_val=0.1, phase_val=1.0) |
| 230 | + model_2d.offset.fixed = False |
| 231 | + model_2d.phase.fixed = False |
| 232 | + |
| 233 | + fitter = Fitter(model_2d, model_2d) |
| 234 | + result = fitter.fit(x=x_2d, y=z_data.ravel()) |
| 235 | +
|
| 236 | +Accessing Fit Results |
| 237 | +-------------------- |
| 238 | + |
| 239 | +The ``FitResults`` object contains comprehensive information about the fit: |
| 240 | + |
| 241 | +.. code-block:: python |
| 242 | +
|
| 243 | + result = fitter.fit(x=x_data, y=y_data) |
| 244 | + |
| 245 | + # Fit statistics |
| 246 | + print(f"Chi-squared: {result.chi2}") |
| 247 | + print(f"Reduced chi-squared: {result.reduced_chi}") |
| 248 | + print(f"Number of parameters: {result.n_pars}") |
| 249 | + print(f"Success: {result.success}") |
| 250 | + |
| 251 | + # Parameter values and uncertainties |
| 252 | + for param_name, value in result.p.items(): |
| 253 | + error = result.errors.get(param_name, 0.0) |
| 254 | + print(f"{param_name}: {value} ± {error}") |
| 255 | + |
| 256 | + # Calculated values and residuals |
| 257 | + y_calculated = result.y_calc |
| 258 | + residuals = result.residual |
| 259 | + |
| 260 | + # Plot results |
| 261 | + import matplotlib.pyplot as plt |
| 262 | + plt.figure(figsize=(10, 4)) |
| 263 | + plt.subplot(121) |
| 264 | + plt.plot(x_data, y_data, 'o', label='Data') |
| 265 | + plt.plot(x_data, y_calculated, '-', label='Fit') |
| 266 | + plt.legend() |
| 267 | + plt.subplot(122) |
| 268 | + plt.plot(x_data, residuals, 'o') |
| 269 | + plt.axhline(0, color='k', linestyle='--') |
| 270 | + plt.ylabel('Residuals') |
| 271 | +
|
| 272 | +Developer Guidelines |
| 273 | +------------------- |
| 274 | + |
| 275 | +Creating Custom Models |
| 276 | +~~~~~~~~~~~~~~~~~~~~~~ |
| 277 | + |
| 278 | +For developers building scientific components: |
| 279 | + |
| 280 | +.. code-block:: python |
| 281 | +
|
| 282 | + from easyscience import ObjBase, Parameter |
| 283 | + |
| 284 | + class CustomModel(ObjBase): |
| 285 | + def __init__(self, param1_val=1.0, param2_val=0.0): |
| 286 | + # Always create Parameters with appropriate bounds and units |
| 287 | + param1 = Parameter("param1", param1_val, min=-10, max=10, unit="m/s") |
| 288 | + param2 = Parameter("param2", param2_val, min=0, max=1, fixed=True) |
| 289 | + |
| 290 | + # Call parent constructor with named parameters |
| 291 | + super().__init__("custom_model", param1=param1, param2=param2) |
| 292 | + |
| 293 | + def __call__(self, x): |
| 294 | + # Implement your model calculation |
| 295 | + return self.param1.value * x + self.param2.value |
| 296 | + |
| 297 | + def get_fit_parameters(self): |
| 298 | + # This is automatically implemented by ObjBase |
| 299 | + # Returns only non-fixed parameters |
| 300 | + return super().get_fit_parameters() |
| 301 | +
|
| 302 | +Best Practices |
| 303 | +~~~~~~~~~~~~~ |
| 304 | + |
| 305 | +1. **Always set appropriate bounds** on parameters to constrain the search space |
| 306 | +2. **Use meaningful units** for physical parameters |
| 307 | +3. **Fix parameters** that shouldn't be optimized |
| 308 | +4. **Test with different minimizers** for robustness |
| 309 | +5. **Validate results** by checking chi-squared and residuals |
| 310 | + |
| 311 | +Error Handling |
| 312 | +~~~~~~~~~~~~~ |
| 313 | + |
| 314 | +.. code-block:: python |
| 315 | +
|
| 316 | + from easyscience.fitting.minimizers import FitError |
| 317 | + |
| 318 | + try: |
| 319 | + result = fitter.fit(x=x_data, y=y_data) |
| 320 | + if not result.success: |
| 321 | + print(f"Fit failed: {result.message}") |
| 322 | + except FitError as e: |
| 323 | + print(f"Fitting error: {e}") |
| 324 | + except Exception as e: |
| 325 | + print(f"Unexpected error: {e}") |
| 326 | +
|
| 327 | +Testing Patterns |
| 328 | +~~~~~~~~~~~~~~~ |
| 329 | + |
| 330 | +When writing tests for fitting code: |
| 331 | + |
| 332 | +.. code-block:: python |
| 333 | +
|
| 334 | + import pytest |
| 335 | + from easyscience import global_object |
| 336 | + |
| 337 | + @pytest.fixture |
| 338 | + def clear_global_map(): |
| 339 | + """Clear global map before each test""" |
| 340 | + global_object.map._clear() |
| 341 | + yield |
| 342 | + global_object.map._clear() |
| 343 | + |
| 344 | + def test_model_fitting(clear_global_map): |
| 345 | + # Create model and test fitting |
| 346 | + model = CustomModel() |
| 347 | + model.param1.fixed = False |
| 348 | + |
| 349 | + # Generate test data |
| 350 | + x_test = np.linspace(0, 10, 50) |
| 351 | + y_test = 2.5 * x_test + 0.1 * np.random.normal(size=len(x_test)) |
| 352 | + |
| 353 | + # Fit and verify |
| 354 | + fitter = Fitter(model, model) |
| 355 | + result = fitter.fit(x=x_test, y=y_test) |
| 356 | + |
| 357 | + assert result.success |
| 358 | + assert model.param1.value == pytest.approx(2.5, abs=0.1) |
| 359 | +
|
| 360 | +This comprehensive guide covers the essential aspects of fitting in EasyScience, from basic usage to advanced developer patterns. |
| 361 | +The examples are drawn from the actual test suite and demonstrate real-world usage patterns. |
| 362 | + |
0 commit comments