lactationcurve

A package for fitting dairy animal lactation curves, evaluating lactation curve characteristics (LCCs) (time to peak, peak yield, cumulative yield, persistency), and computing 305-day milk yield using the ICAR guideline.

Contact: Meike van Leerdam, mbv32@cornell.edu

Authors: Meike van Leerdam, Douwe de Kok, Judith Osei-Tete, Lucia Trapanese

Initial authored: 2025‑08‑12

Updated: 2026-06-03


The 305 day yield for milk, fat, and protein is a widely used metric in dairy production, and the International Committee for Animal Recording (ICAR) provides guidelines outlining approved methods for its calculation. However, a global survey of milk recording organizations revealed substantial variation in how these methods are implemented. The Test Interval Method is used by 74% of the organizations, reflecting a preference for methodological simplicity, but it comes with trade-offs in estimation accuracy. The use of the other approved methods showed wide variation in correction factors, standard lactation curves, test day definitions, minimum sample requirements, and exclusion criteria. Such inconsistencies can introduce yield variability that complicates comparisons, for example in international breeding value evaluation, and limit the metric’s usefulness in universal models, such as decision support tools. Thus, the objective of this work was to reformulate the ICAR guideline section 2, procedure 2, into a unified, transparent, and accessible software implementation to improve standardization, enhance documentation, support continuous development, and increase the accuracy of 305 day yield estimation.

To achieve this, the ICAR guideline was translated into an open-source Python package that serves as a reference implementation for 305-day yield calculation, with lactation curve modelling at its core. In addition to the methods described in the original guideline, the package incorporates 14 lactation curve models, including both traditional and Bayesian fitting approaches as well as AI-based models. It also provides tools to derive biologically relevant characteristics such as time to peak, peak yield, cumulative yield, and persistency.

The framework can be directly integrated into analytical workflows, allowing users to calculate 305-day yields, fit and compare lactation curves, and derive key lactation characteristics through a single function call. These functionalities are further supported through an interactive website and an openly available GitHub repository. In addition, the project includes an online validation platform that enables users to use standardized lactation data to compare self-estimated 305-day yields against both reference calculations (from the package) and observed daily milk yields.

We encourage everyone to use, test, and contribute to the package, which is available under the AGPL-3.0-or-later license. We welcome feedback and suggestions for improvement, and we are committed to maintaining and updating the package to ensure it remains a valuable resource for the dairy industry and research community. For bug reports or feature requests, please submit them on the GitHub issues page or contact us directly via email.

Main Lactation curve models implemented:

Fischer: A combination of an exponential and a linear component, characterized by an early exponential phase followed by an approximately linear decline.

Wood: An incomplete gamma-type function combining a power-law growth term t^b with an exponential decay term e^(-ct).

Ali & Schaeffer: A linear regression model based on quadratic polynomials in standardized time (t/305) and in the log-transformed inverse time log⁡(305/t), where the log term increases flexibility in modelling the ascending phase and peak of lactation.

Wilmink: A combination of an exponential and a linear component. In contrast to the Fischer model, the additional parameter scaling the exponential term increases flexibility in describing early-lactation dynamics and peak formation.

MilkBot: An empirical, mechanistically motivated four-parameter model describing lactation as the development and decay of udder capacity. It consists of a ramp-up phase and exponential decline, with parameters controlling scale, onset, growth rate, and decay. Both Bayesian and frequentist fitting approaches are implemented for the MilkBot model, allowing users to choose between traditional optimization methods and a probabilistic framework that incorporates prior knowledge.

Additional models available for a.o. symbolic lactation curve characteristics (LCC) derivations: Brody, Sikka, Nelder, Dhanoa, Emmans, Hayashi, Rook, Dijkstra, Prasad.


Model Formulas

  • Wood : y(t) = a * t^b * exp(-c * t)
  • Wilmink : y(t) = a + b * t + c * exp(k * t) with default k = -0.05
  • Ali & Schaeffer : t_scaled = t / 305, L = ln(305 / t) y(t) = a + b*t_scaled + c*t_scaled^2 + d*L + k*L^2
  • Fischer : y(t) = a - b*t - a*exp(-c*t)
  • MilkBot : y(t) = a * (1 - exp((c - t)/b) / 2) * exp(-d*t)

Features

  • Frequentist fitting (numeric optimization):
    • Wood, Wilmink, Ali & Schaeffer, Fischer, MilkBot
  • Frequentist fitting (algebraic least squares):
    • MilkBot
  • Bayesian fitting via MilkBot API:
    • MilkBot
  • Lactation Curve Characteristics — symbolic + numeric:
    • time_to_peak, peak_yield, cumulative_milk_yield, persistency
  • ICAR procedures cumulative milk yield:
    • Test Interval Method
    • Interpolation Standard Lactation Curve (ISLC) Method
    • Best Predict Method
  • Input validation/normalization via validate_and_prepare_inputs
  • Caching of symbolic expressions for performance

API Overview

The package is organized into three main modules:

  1. lactationcurve.fitting
  2. lactationcurve.characteristics
  3. lactationcurve.preprocessing

Output Types Summary Of Most Important Functions

| Function | Output |

|---------|--------|

| fit_lactation_curve | Predicted yields (np.ndarray) |

| get_lc_parameters | Tuple of numerical parameters |

| bayesian_fit_milkbot_single_lactation | Dict of MilkBot parameters |

| lactation_curve_characteristic_function | (expr, params, func) |

| calculate_characteristic | float (LCC value) |

| test_interval_method | DataFrame with 305‑day totals per TestId |

| ISLC_method | DataFrame with 305‑day totals per TestId |

| best_predict_method | DataFrame with 305‑day totals per TestId |


The meaning of a TestId

The TestId is an identifier for a lactation, which can be used to group records belonging to the same lactation together. It is not the same as a cow ID, as a cow can have multiple lactations (e.g., across different calvings). If a TestId column is not provided, the package will assume all records belong to a single lactation and will create a TestId column with all values set to 0.


Often used abreviations

  • DIM: Days in Milk (the days since calving in the current lactation)
  • LC: Lactation Curve
  • LCC: Lactation Curve Characteristic
  • ISLC: Interpolation using the Standard Lactation Curve
  • ICAR: International Committee for Animal Recording
  • API: Application Programming Interface

Bayesian Fitting (MilkBot API)

  • Set fitting="bayesian" and model="milkbot" in fit_lactation_curve or calculate_characteristic.
  • Provide an API key via .env
  • Choose priors via custom_priors:

    • "CHEN" → Chen et al. 2023 published priors
    • dict → Custom priors in MilkBot format (overrides continent) Custom priors have a specific format, to help you build them, use the build_prior helper function. To make your own prior you need for each MilkBot parameter (scale,ramp,decay,onset) to specify a mean and a standard deviation (std). Also provide a standard deviation for milk yield through seMilk to specify the expected noise in the data. This seMilk is default set to 4 kg to reflect typical day-to-day variation in milk yield, but you can adjust it based on the expected variability in your data.
  • if no custom priors are provided, the default MilkBot priors will be used: In that case set cutsom_priors to None and specify the desired continent and breed. continent options:

    • "USA" → MilkBot USA priors (default)
    • "EU" → MilkBot EU priors > mainly estimates lower milk production

    breed options:

    • H → Holstein (default)
    • J → Jersey

    If also parity is provided, the continent-specific priors will be further refined by parity-specific priors. The priors are sensitive to the used metric. The default is kg, but if you use lb, specify milk_unit="lb" to use the appropriate priors.

  • The helper bayesian_fit_milkbot_single_lactation(...) normalizes differing API responses.

  • The key can be requested by sending an email to Jim Ehrlich jehrlich@MilkBot.com.
  • More information about the API can be found in the API documentation, or in the corresponding paper.


Tutorials

Hands-on programming example notebooks are available through GitHub tutorials and example notebooks. The tutorials will guide you through how to work with the package, including how to import the package, fit lactation curves, calculate 305-day yields, and derive lactation curve characteristics.


Citing the lactationcurve package

If you use the lactationcurve package in your research, please consider citing it as follows:

van Leerdam, M. B., de Kok, D., Osei-Tete, J. A., & Hostens, M. (2026). Bovi-analytics/bovi: v.1.1.6. (v.1.1.6). Zenodo. https://doi.org/10.5281/zenodo.18715145

If you also use the Bayesian fitting functionality that relies on the MilkBot API, please also cite the following paper:

Ehrlich, J.L., 2013. Quantifying inter-group variability in lactation curve shape and magnitude with the MilkBot lactation model. PeerJ 1, e54. https://doi.org/10.7717/peerj.54

If you use the 305-day yield calculation methods based on the ICAR guideline, please also cite the following paper: Best Predict method: VanRaden, P. M. (1997). Lactation yields and accuracies computed from test day yields and (co) variances by best prediction. Journal of dairy science, 80(11), 3015-3022.

ISLC: Wilmink, J. B. M. (1987). Comparison of different methods of predicting 305-day milk yield using means calculated from within-herd lactation curves. Livestock Production Science, 17, 1-17.

Test Interval Method: Sargent, F. D., V. H. Lyton, and 0. G. Wall, J r . 1968. Test interval method of calculating Dairy Herd Improvement Association records. Journal of dairy science, 51-170.


License

GNU AGPLv3 License


Current version

0.1.dev1+g7b9c117e1

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