lactationcurve

  1"""
  2
  3A package for fitting **dairy animal lactation curves**, evaluating
  4**lactation curve characteristics (LCCs)** (time to peak, peak yield,
  5cumulative yield, persistency), and computing **305-day milk yield**
  6using the **ICAR guideline**.
  7
  8> **Contact:** Meike van Leerdam, mbv32@cornell.edu
  9>
 10> **Authors:** Meike van Leerdam, Douwe de Kok, Judith Osei-Tete, Lucia Trapanese
 11
 12> **Initial authored:** 2025‑08‑12
 13
 14> **Updated:** 2026-06-03
 15
 16---
 17
 18  The 305 day yield for milk, fat, and protein is a widely used metric
 19  in dairy production, and the International Committee for Animal
 20  Recording (ICAR) provides guidelines outlining approved methods for
 21  its calculation. However, a global survey of milk recording
 22  organizations revealed substantial variation in how these methods are
 23  implemented. The Test Interval Method is used by 74% of the
 24  organizations, reflecting a preference for methodological simplicity,
 25  but it comes with trade-offs in estimation accuracy. The use of the
 26  other approved methods showed wide variation in correction factors,
 27  standard lactation curves, test day definitions, minimum sample
 28  requirements, and exclusion criteria. Such inconsistencies can
 29  introduce yield variability that complicates comparisons, for
 30  example in international breeding value evaluation, and limit the
 31  metric’s usefulness in universal models, such as decision support
 32  tools. Thus, the objective of this work was to reformulate the ICAR
 33  guideline section 2, procedure 2, into a unified, transparent, and
 34  accessible software implementation to improve standardization,
 35  enhance documentation, support continuous development, and increase
 36  the accuracy of 305 day yield estimation.
 37
 38  To achieve this, the ICAR guideline was translated into an
 39  open-source Python package that serves as a reference implementation
 40  for 305-day yield calculation, with lactation curve modelling at its
 41  core. In addition to the methods described in the original guideline,
 42  the package incorporates 14 lactation curve models, including both
 43  traditional and Bayesian fitting approaches as well as AI-based
 44  models. It also provides tools to derive biologically relevant
 45  characteristics such as time to peak, peak yield, cumulative yield,
 46  and persistency.
 47
 48  The framework can be directly integrated into analytical workflows,
 49  allowing users to calculate 305-day yields, fit and compare
 50  lactation curves, and derive key lactation characteristics through a
 51  single function call. These functionalities are further supported
 52  through an interactive website and an openly available [GitHub
 53  repository](https://github.com/Bovi-analytics/bovi).
 54  In addition, the project includes an online validation
 55  platform that enables users to use standardized lactation data to
 56  compare self-estimated 305-day yields against both reference
 57  calculations (from the package) and observed daily milk yields.
 58
 59  We encourage everyone to use, test, and contribute to the package,
 60  which is available under the MIT license. We welcome feedback and
 61  suggestions for improvement, and we are committed to maintaining and
 62  updating the package to ensure it remains a valuable resource
 63  for the dairy industry and research community.
 64  For bug reports or feature requests, please submit them on the
 65  [GitHub issues page](https://github.com/Bovi-analytics/bovi/issues)
 66  or contact us directly via [email](mbv32@cornell.edu).
 67
 68
 69## Main Lactation curve models implemented:
 70
 71**Fischer**:
 72A combination of an exponential and a linear component, characterized by an
 73early exponential phase followed by an approximately linear decline.
 74
 75**Wood**:
 76An incomplete gamma-type function combining a power-law growth term t^b with
 77an exponential decay term e^(-ct).
 78
 79**Ali & Schaeffer**:
 80A linear regression model based on quadratic polynomials in standardized time
 81(t/305) and in the log-transformed inverse time log⁡(305/t), where the log
 82term increases flexibility in modelling the ascending phase and peak of
 83lactation.
 84
 85**Wilmink**:
 86A combination of an exponential and a linear component. In contrast to the
 87Fischer model, the additional parameter scaling the exponential term increases
 88flexibility in describing early-lactation dynamics and peak formation.
 89
 90**MilkBot**:
 91An empirical, mechanistically motivated four-parameter model describing
 92lactation as the development and decay of udder capacity. It consists of a
 93ramp-up phase and exponential decline, with parameters controlling scale,
 94onset, growth rate, and decay.
 95Both Bayesian and frequentist fitting approaches are implemented for the
 96MilkBot model, allowing users to choose between traditional optimization
 97methods and a probabilistic framework that incorporates prior knowledge.
 98
 99Additional models available for a.o. symbolic lactation curve
100characteristics (LCC) derivations:
101**Brody**, **Sikka**, **Nelder**, **Dhanoa**, **Emmans**,
102**Hayashi**, **Rook**, **Dijkstra**, **Prasad**.
103
104---
105
106## Model Formulas
107
108* **Wood** : `y(t) = a * t^b * exp(-c * t)`
109* **Wilmink** : `y(t) = a + b * t + c * exp(k * t)` with default `k = -0.05`
110* **Ali & Schaeffer** :  `t_scaled = t / 305`, `L = ln(305 / t)`
111  `y(t) = a + b*t_scaled + c*t_scaled^2 + d*L + k*L^2`
112* **Fischer** : `y(t) = a - b*t - a*exp(-c*t)`
113* **MilkBot** : `y(t) = a * (1 - exp((c - t)/b) / 2) * exp(-d*t)`
114
115---
116
117## Features
118
119- **Frequentist fitting** (numeric optimization):
120  - Wood, Wilmink, Ali & Schaeffer, Fischer, MilkBot
121- **Frequentist fitting** (algebraic least squares):
122  - MilkBot
123- **Bayesian fitting via MilkBot API**:
124  - MilkBot
125- **Lactation Curve Characteristics** — symbolic + numeric:
126  - time_to_peak, peak_yield, cumulative_milk_yield, persistency
127- **ICAR procedures cumulative milk yield:**
128  - Test Interval Method
129  - Interpolation Standard Lactation Curve (ISLC) Method
130  - Best Predict Method
131- Input validation/normalization via `validate_and_prepare_inputs`
132- Caching of symbolic expressions for performance
133
134---
135
136## API Overview
137
138The package is organized into three main modules:
139
1401. `lactationcurve.fitting`
1412. `lactationcurve.characteristics`
1423. `lactationcurve.preprocessing`
143
144---
145
146## Output Types Summary Of Most Important Functions
147
148| Function | Output |
149
150|---------|--------|
151
152| `fit_lactation_curve` | Predicted yields (np.ndarray) |
153
154| `get_lc_parameters` | Tuple of numerical parameters |
155
156| `bayesian_fit_milkbot_single_lactation` | Dict of MilkBot parameters |
157
158| `lactation_curve_characteristic_function` | (expr, params, func) |
159
160| `calculate_characteristic` | float (LCC value) |
161
162| `test_interval_method` | DataFrame with 305‑day totals per TestId |
163
164| `ISLC_method` | DataFrame with 305‑day totals per TestId  |
165
166| `best_predict_method` | DataFrame with 305‑day totals per TestId |
167
168---
169
170## The meaning of a TestId
171
172
173The `TestId` is an identifier for a lactation,
174which can be used to group records belonging to the same lactation together.
175It is not the same as a cow ID, as a cow can have multiple lactations
176(e.g., across different calvings).
177If a `TestId` column is not provided,
178the package will assume all records belong to a single lactation
179and will create a `TestId` column with all values set to 0.
180
181---
182
183## Often used abreviations
184
185- **DIM**: Days in Milk (the days since calving in the current lactation)
186- **LC**: Lactation Curve
187- **LCC**: Lactation Curve Characteristic
188- **ISLC**: Interpolation using the Standard Lactation Curve
189- **ICAR**: International Committee for Animal Recording
190- **API**: Application Programming Interface
191
192
193---
194
195## Bayesian Fitting (MilkBot API)
196
197* Set `fitting="bayesian"` and `model="milkbot"` in
198  `fit_lactation_curve` or `calculate_characteristic`.
199* Provide an **API key** via .env
200* Choose priors via custom_priors:
201    - "[CHEN](https://github.com/Bovi-analytics/Chen-et-al-2023b)" → Chen et al. 2023
202      published priors
203    - dict    → Custom priors in MilkBot format (overrides `continent`)
204      Custom priors have a specific format, to help you build them,
205      use the `build_prior` helper function. To make your own prior you need
206      for each MilkBot parameter
207      (scale,ramp,decay,onset) to specify a mean and a standard deviation (std).
208      Also provide a standard deviation for milk yield through seMilk to
209      specify the expected noise in the data.
210      This seMilk is default set to 4 kg to reflect typical day-to-day variation in milk yield,
211      but you can adjust it based on the expected variability in your data.
212
213* if no custom priors are provided, the default MilkBot priors will be used:
214    In that case set cutsom_priors to None and specify the desired continent and breed.
215    continent options:
216    - "USA"   → MilkBot USA priors (default)
217    - "EU"    → MilkBot EU priors > mainly estimates lower milk production
218
219    breed options:
220    - H → Holstein (default)
221    - J → Jersey
222
223    If also parity is provided, the continent-specific priors will be
224    further refined by parity-specific priors.
225    The priors are sensitive to the used metric. The default is kg, but if you use lb,
226    specify milk_unit="lb" to use the appropriate priors.
227
228* The helper `bayesian_fit_milkbot_single_lactation(...)`
229  normalizes differing API responses.
230* The key can be requested by sending an email to Jim Ehrlich
231  [jehrlich@MilkBot.com](mailto:jehrlich@MilkBot.com).
232* More information about the API can be found in the
233  [API documentation](https://api.milkbot.com/), or in the
234  corresponding
235  [paper](https://peerj.com/articles/54/#MainContent).
236
237---
238
239## Citing the lactationcurve package
240
241If you use the `lactationcurve` package in your research, please consider citing it as follows:
242
243*van Leerdam, M. B., de Kok, D., Osei-Tete, J. A., &
244Hostens, M. (2026). Bovi-analytics/bovi:
245v.1.1.6. (v.1.1.6). Zenodo.
246https://doi.org/10.5281/zenodo.18715145*
247
248
249If you also use the Bayesian fitting functionality that relies
250on the MilkBot API, please also cite the following paper:
251
252*Ehrlich, J.L., 2013. Quantifying inter-group variability
253in lactation curve shape and magnitude with the MilkBot
254lactation model. PeerJ 1, e54.
255https://doi.org/10.7717/peerj.54*
256
257If you use the 305-day yield calculation methods based on the ICAR guideline,
258please also cite the following paper:
259Best Predict method:
260*VanRaden, P. M. (1997). Lactation yields and accuracies computed from test
261day yields and (co) variances by best prediction.
262Journal of dairy science, 80(11), 3015-3022.*
263
264ISLC:
265*Wilmink, J. B. M. (1987).
266Comparison of different methods of predicting 305-day milk yield using means
267calculated from within-herd lactation curves. Livestock Production Science, 17, 1-17.*
268
269Test Interval Method:
270*Sargent, F. D., V. H. Lyton, and 0. G. Wall, J r . 1968.
271Test interval method of calculating Dairy Herd Improvement Association records.
272Journal of dairy science, 51-170.*
273
274---
275
276## License
277
278[MIT License](https://github.com/Bovi-analytics/lactation_curve_core/blob/master/LICENSE)
279
280
281---
282
283## Version v.1.1.6
284
285"""
286
287
288# import submodules to make them available at the package level
289
290from . import characteristics, fitting, preprocessing
291
292__all__ = ["fitting", "characteristics", "preprocessing"]
293# from .characteristics import (
294#     calculate_characteristic,
295#     lactation_curve_characteristic_function,
296#     numeric_cumulative_yield,
297#     numeric_peak_yield,
298#     numeric_time_to_peak,
299#     persistency_fitted_curve,
300#     persistency_milkbot,
301#     persistency_wood,
302#     test_interval_method,
303# )
304# from .fitting import (
305#     ali_schaeffer_model,
306#     bayesian_fit_milkbot_single_lactation,
307#     brody_model,
308#     dhanoa_model,
309#     dijkstra_model,
310#     emmans_model,
311#     fischer_model,
312#     fit_lactation_curve,
313#     get_chen_priors,
314#     get_lc_parameters,
315#     get_lc_parameters_least_squares,
316#     hayashi_model,
317#     milkbot_model,
318#     nelder_model,
319#     prasad_model,
320#     rook_model,
321#     sikka_model,
322#     wilmink_model,
323#     wood_model,
324#     build_prior,
325# )
326# from .preprocessing import (
327#     PreparedInputs,
328#     standardize_lactation_columns,
329#     validate_and_prepare_inputs,
330# )
331
332# __all__ = [
333#     # Preprocessing
334#     "PreparedInputs",
335#     "standardize_lactation_columns",
336#     "validate_and_prepare_inputs",
337#     # Fitting
338#     "ali_schaeffer_model",
339#     "bayesian_fit_milkbot_single_lactation",
340#     "brody_model",
341#     "dhanoa_model",
342#     "dijkstra_model",
343#     "emmans_model",
344#     "fischer_model",
345#     "fit_lactation_curve",
346#     "get_chen_priors",
347#     "get_lc_parameters",
348#     "get_lc_parameters_least_squares",
349#     "hayashi_model",
350#     "milkbot_model",
351#     "nelder_model",
352#     "prasad_model",
353#     "rook_model",
354#     "sikka_model",
355#     "wilmink_model",
356#     "wood_model",
357#     # Characteristics
358#     "calculate_characteristic",
359#     "lactation_curve_characteristic_function",
360#     "numeric_cumulative_yield",
361#     "numeric_peak_yield",
362#     "numeric_time_to_peak",
363#     "persistency_fitted_curve",
364#     "persistency_milkbot",
365#     "persistency_wood",
366#     "test_interval_method",
367#     "build_prior",
368# ]
369
370# Expose package version (try metadata, fall back to a sensible dev string)
371try:
372    from importlib.metadata import PackageNotFoundError, version
373except Exception:
374    try:
375        from importlib_metadata import PackageNotFoundError, version  # type: ignore
376    except Exception:
377        version = None
378        PackageNotFoundError = Exception
379
380if version:
381    try:
382        __version__ = version("lactationcurve")
383    except PackageNotFoundError:
384        __version__ = "0+dev"
385else:
386    __version__ = "0+dev"
387
388__all__.append("__version__")