The phenology module

The phenology module.

The FAPAR limitation module

Class to compute the fAPAR_max and annual peak Leaf Area Index (LAI).

The FaparLimitation class and the FaparLimitation.from_pmodel() are designed to work with inputs that can have multiple dimensions. The first axis is _always_ assumed to represent a time series of annual observations. If the inputs are one dimensional, then this is a time series for a single site; if they are three dimensional then these are observations for a grid of sites. Usually all array inputs will have the same shape but note the following instances where you might need to take care with array broadcasting.

  • Growing season length might well be constant across sites. If so - for example

    with 3D data - the input would need shape (N, 1, 1) to broadcast N years of data over the array of sites.

  • The climatological aridity index is very likely to be constant through time.

    If so, the array should have a singleton first dimension to broadcast an observation per site across observations. For example, with (10, 3, 3) data (10 years for a 3x3 grid of sites), the aridity index could be provided as (1,3,3) to broadcast the aridity index across each year. It could also use a single scalar value to use the same aridity index for all sites.

class FaparLimitation(annual_total_potential_gpp: ndarray[tuple[Any, ...], dtype[floating]], annual_mean_ca: ndarray[tuple[Any, ...], dtype[floating]], annual_mean_chi: ndarray[tuple[Any, ...], dtype[floating]], annual_mean_vpd: ndarray[tuple[Any, ...], dtype[floating]], annual_total_precip: ndarray[tuple[Any, ...], dtype[floating]], annual_growing_season_length: ndarray[tuple[Any, ...], dtype[floating]], years: ndarray[tuple[Any, ...], dtype[datetime64]], method: str = 'cai', phenology_const: PhenologyConstNew = PhenologyConstNew(), **kwargs: ndarray[tuple[Any, ...], dtype[floating]])

Compute maximum annual fAPAR and LAI.

Experimental

Be aware that FaparLimitation is an experimental feature and the API and any calculated values may change between major releases.

This class calculates maximum annual fAPAR, which can be limited either by the availability of light energy ($f_{APAR_{c}}$) or by the availability of water ($f_{APAR_{w}}$). The equations for these two variables, following (Cai et al., 2025) are:

\[ \begin{align*} f_{APAR_{c}} &= 1 - \frac{z}{k A_0}\\ f_{APAR_{w}} &= \left(\frac{ c_a \left( 1 - \chi \right)}{ 1.6 D }\right) \left(\frac{ f_0 P }{ A_0 }\right) \\ \end{align*} \]

The maximum fAPAR is then calculated as a function of \(f_{APAR_{c}}\) and \(f_{APAR_{w}}\). In these equations:

  • \(z\) accounts for the growth and maintenance costs of leaves.

  • \(f_0\) accounts for water limitation on annual assimulation and is is the ratio of annual total transpiration to annual total precipitation.

The other variables are the required arguments to the class defined below.

There are different approaches to estimating \(z\) and \(f_0\) and to calculating the maximum fAPAR from the two inputs. For details see:

The maximum annual LAI can then be calculated using Beer’s law as:

\[L_{max} = - ( 1 / k ) \ln {1 -f_{APAR_{max}}}\]

The most common source of most of the variables needed to calculate maximum fAPAR is a P Model, and the from_pmodel() method can be used to estimate maximum fAPAR directly from a fitted P Model.

Parameters:

  • annual_total_potential_gpp – The annual sum of potential GPP (\(A_0, \text{mol C m}^{-2} \text{year}^{-1}\))

  • annual_mean_ca – The ambient CO2 partial pressure during the growing season (\(c_a\), Pa)

  • annual_mean_chi – The annual mean ratio of ambient to leaf CO2 partial during the growing season (\(\chi\), Pa)

  • annual_mean_vpd – The annual mean vapour pressure deficit during the growing season (\(D\), Pa)

  • annual_total_precip – The annual total precipitation (\(P, \text{mol m}^{-2} \text{year}^{-1}\))

  • annual_growing_season_length – The length of the growing season in days for each year (\(G\), days)

  • years – An array of year datetimes for the observations.

  • method – The method to be applied when calculating maximum fAPAR, defaulting to

  • cai.

  • phenology_const – An instance of PhenologyConstNew

  • **kwargs – Any additional variables required by specific method choices.

classmethod from_pmodel(pmodel: PModelABC, growing_season: ndarray[tuple[Any, ...], dtype[bool]], precip: ndarray[tuple[Any, ...], dtype[floating]], datetimes: ndarray[tuple[Any, ...], dtype[datetime64]] | None = None, method: str = 'cai', phenology_const: PhenologyConstNew = PhenologyConstNew(), **kwargs: ndarray[tuple[Any, ...], dtype[floating]]) → FaparLimitation

Create a FaparLimitation instance from a P Model and other inputs.

The annual summary values of \(A_0, c_a, \chi\) and \(D\) used by the FaparLimitation() class can be taken directly from the predictions of a P Model. This method automatically extracts the required data from a fitted P Model and returns a FaparLimitation instance.

Note

The calculation of fAPAR limitation requires estimates of potential GPP, so the PModelEnvironment instance used to fit the model must set fapar to be one.

Some additional information is needed:

  • The calculation requires annual summaries of variables, so the datetimes argument must be used to provide an array of datetimes for each observation.

  • The annual mean values \(c_a, \chi\) and \(D\) should be estimated during the growing season, so the growing_season argument must be used to provide an array of boolean values indicating which observations should be treated as in the growing season.

  • The calculation requires estimates of precipitation, so the precipitation argument must provide estimates of total precipitation during each observations in moles of water per metre squared.

  • The calculation of the \(f_0\) parameter requires estimates of site specific aridity index.

The method accepts both standard and subdaily P Models and automatically uses the actual time intervals between observations to calculate the required weighted annual means and sums. This might lead to unexpected values: the yearly mean of monthly values \(1, 2, \dots, 12\) would not be 6.5 because the monthly values are weighted according to the length of the month.

Parameters:

  • pmodel – A pyrealm.pmodel.pmodel.PModel or pyrealm.pmodel.pmodel.SubdailyPModel instance, fitted with fapar fixed at one.

  • datetimes – An array giving the datetimes of observations.

  • growing_season – A boolean array indicating which observations are to be considered as part of the growing season.

  • precip – An array of precipitation for each observation.

  • method – The method to be used in calculating maximum fAPAR, defaulting to cai.

  • phenology_const – An instance of PhenologyConstNew

  • **kwargs – Any additional variables required by specific method choices.

summarize(dp: int = 2) → None

Print summary of estimates of fAPAR limitation.

Prints a summary of the calculated values in a FaparLimitation instance including the mean, range and number of nan values.

Parameters:

dp – The number of decimal places used in rounding summary stats.

annual_growing_season_length

Annual growing season length (\(G\), days)

annual_mean_ca

Ambient CO2 partial pressure during the growing season (\(c_a\), Pa)

annual_mean_chi

Annual mean ratio of ambient to leaf CO2 partial during the growing season (\(\chi\), Pa)

annual_mean_vpd

Annual mean vapour pressure deficit during the growing season (\(D\), Pa)

annual_total_potential_gpp

The annual sum of potential GPP (\(A_0, \text{mol C m}^{-2} \text{year}^{-1}\))

annual_total_precip

Annual total precipitation (\(P, \text{mol m}^{-2} \text{year}^{-1}\))

lai_max: ndarray[tuple[Any, ...], dtype[floating]]

Estimated annual maximum LAI (unitless)

years

The year of each observation.

class FaparLimitationMethodABC(fapar_limitation: FaparLimitation, phenology_const: PhenologyConstNew)

Abstract base class for methods computing maximum annual fAPAR and LAI.

Experimental

Be aware that FaparLimitationMethodABC is an experimental feature and the API and any calculated values may change between major releases.

abstract calculate_maximum_fapar(energy_limited_fapar: ndarray[tuple[Any, ...], dtype[_ScalarT]], water_limited_fapar: ndarray[tuple[Any, ...], dtype[_ScalarT]]) → ndarray[tuple[Any, ...], dtype[floating]]

Calculate the maximum fAPAR.

Provides the method specific calculation of maximum fAPAR from the energy and water limited maximum values and should return the calculated maximum fAPAR.

Parameters:

  • energy_limited_fapar – The maximum fAPAR given energy limitation.

  • water_limited_fapar – The maximum fAPAR given water limitation.

abstract set_z_and_f0() → None

Sets the f0 and z values for the method.

attrs: tuple[tuple[str, str], ...]

A tuple of attributes to be reported in the FaparLimitation.summarize() output for the method.

method: str

A short method name used to identify the class in FAPAR_LIMITATION_METHOD_CLASS_REGISTRY.

requires: tuple[str, ...]

A tuple of any additional variables required to create a FaparLimitation instance using the method.

class FaparLimitationMethodCai(fapar_limitation: FaparLimitation, phenology_const: PhenologyConstNew)

Compute maximum annual fAPAR and LAI following Cai et al. (2025).

Experimental

Be aware that FaparLimitationMethodCai is an experimental feature and the API and any calculated values may change between major releases.

The method of Cai et al. (2025) uses a single global value of \(z\) but models \(f_0\) as a function of site-specific aridity, expressed as the climatological ratio PET/P, (see set_z_and_f0()). The annual maximum fAPAR is then the simple minimum of the site values for water and energy limited fAPAR.

calculate_maximum_fapar(energy_limited_fapar: ndarray[tuple[Any, ...], dtype[_ScalarT]], water_limited_fapar: ndarray[tuple[Any, ...], dtype[_ScalarT]]) → ndarray[tuple[Any, ...], dtype[_ScalarT]]

Calculate the maximum fAPAR.

The Cai method uses the simple minimum of the water and energy limited fAPAR values as the maximum possible fAPAR.

Parameters:

  • energy_limited_fapar – The maximum fAPAR given energy limitation.

  • water_limited_fapar – The maximum fAPAR given water limitation.

set_z_and_f0() → None

Set the \(z\) and \(f_0\) parameters.

The value \(f_0\) is the ratio of annual total transpiration of annual total precipitation. It is calculated from site specific estimates of the climatological aridity index, calculated as the long term (typically 20 years) total PET over total precipitation (\(AI\), unitless) as:

\[f_0 = a \exp{\left(-b \left(\frac{AI}{c}\right)^2\right)}\]

where \(a,b,c\) are defined in the cai_f0_coefficients attribute.

attrs: tuple[tuple[str, str], ...] = (('lai_max', '-'), ('fapar_max', '-'))

A tuple of attributes to be reported in the FaparLimitation.summarize() output for the method.

energy_limited: ndarray[tuple[Any, ...], dtype[bool]]

Boolean array showing if annual \(fAPAR_{max}\) is water or energy limited.

method: str = 'cai'

A short method name used to identify the class in FAPAR_LIMITATION_METHOD_CLASS_REGISTRY.

requires: tuple[str, ...] = ('aridity_index',)

A tuple of any additional variables required to create a FaparLimitation instance using the method.

class FaparLimitationMethodZhu(fapar_limitation: FaparLimitation, phenology_const: PhenologyConstNew)

Compute maximum annual fAPAR and LAI following Zhu et al. (2026).

Experimental

Be aware that FaparLimitationMethodZhu is an experimental feature and the API and any calculated values may change between major releases.

The method of Zhu et al. (2026) uses single global values of both \(z\) and \(f_0\). However, the annual maximum fAPAR is then a continuous function of water and energy limited fAPAR following a Budyko curve (see calculate_maximum_fapar() for details).

calculate_maximum_fapar(energy_limited_fapar: ndarray[tuple[Any, ...], dtype[_ScalarT]], water_limited_fapar: ndarray[tuple[Any, ...], dtype[_ScalarT]]) → ndarray[tuple[Any, ...], dtype[_ScalarT]]

Calculate the maximum fAPAR.

The Zhu method calculates the maximum fAPAR as a function of the energy and water limited fAPAR, following a function like the Budyko curve (Roderick and Farquhar, 2011):

\[ \begin{align*} c_w &= \frac{f_{APAR_{c}}}{f_{APAR_{w}}} \\ f_{APAR_{max}} = f_{APAR_{w}} \left(1 + c_w - \left( 1 + c_w ^{\bar\omega} \right)^{1/{\bar\omega}} \end{align*} \],

where \(\bar\omega\) is defined in the zhu_budyko constants parameter.

Parameters:

  • energy_limited_fapar – The maximum fAPAR given energy limitation.

  • water_limited_fapar – The maximum fAPAR given water limitation.

set_z_and_f0() → None

Set the \(z\) and \(f_0\) parameters.

This method has fixed values for \(z\) (see zhu_f0) and \(f_0\) (see zhu_z)

attrs: tuple[tuple[str, str], ...] = (('lai_max', '-'), ('fapar_max', '-'))

A tuple of attributes to be reported in the FaparLimitation.summarize() output for the method.

method: str = 'zhu'

A short method name used to identify the class in FAPAR_LIMITATION_METHOD_CLASS_REGISTRY.

requires: tuple[str, ...] = ()

A tuple of any additional variables required to create a FaparLimitation instance using the method.

FAPAR_LIMITATION_METHOD_CLASS_REGISTRY: dict[str, type[FaparLimitationMethodABC]] = {'cai': <class 'pyrealm.phenology.fapar_limitation.FaparLimitationMethodCai'>, 'zhu': <class 'pyrealm.phenology.fapar_limitation.FaparLimitationMethodZhu'>}

A registry for optimal chi calculation classes.

Different implementations of the calculation of optimal chi must all be subclasses of OptimalChiABC abstract base class. This dictionary is used as a registry for defined subclasses and a method name is used to retrieve a particular implementation from this registry. For example:

prentice14_opt_chi = OPTIMAL_CHI_CLASS_REGISTRY['prentice14']

The Phenology module

Computing phenological time series of Leaf Area Index.

This module provides classes to implement different approaches to calculating daily predictions of leaf area index (LAI) from daily time series of potential assimilation. The Phenology class acts as a wrapper around alternative implementations of the calculations, which are implemented as derived subclasses of the PhenologyMethodABC abstract base class.

The methods all require values of annual maximum fAPAR, calculated using FaparLimitation and then a time series of daily assimilation, although different methods can also require additional inputs. The methods then broadly work by calculating a time series of steady state LAI as a function of daily assimilation and then applying a lag, capturing the speed with which plants can put on leaf area, to generate a time series of realised LAI.

Values of daily assimilation are typically going to be taken from a P Model, so the Phenology.from_pmodel() method is provided to calculate daily assimilation from a P Model input and then use this, along with a FaparLimitation instance, to generate an LAI time series.

The Phenology class and the Phenology.from_pmodel() are designed to work with inputs that can have multiple dimensions. The first axis is _always_ assumed to represent a time series of daily observations of potential assimilation. If all the arrays are one dimensional, then this is a time series for a single site; if they are three dimensional then these are observations for a grid of sites.

Usually all array inputs will have the same shape. However, where methods use a climatological measure of aridity, then this is very likely to be constant through time. If so, the array should have a singleton first dimension to broadcast an observation per site across observations. For example, with (365, 3, 3) data (one year of observations for a 3x3 grid of sites), the aridity index could be provided as (1,3,3) to broadcast the aridity index across all observations. It could also use a single scalar value to use the same aridity index for all sites.

class Phenology(fapar_limitation: FaparLimitation, daily_potential_assimilation: ndarray[tuple[Any, ...], dtype[floating]], datetimes: ndarray[tuple[Any, ...], dtype[datetime64]], method: str = 'zhou', **kwargs: Any)

Estimating daily time series of leaf area index (LAI).

Experimental

Be aware that Phenology is an experimental feature and the API and any calculated values may change between major releases.

The maximum fAPAR and LAI for a year (see FaparLimitation) can be combined with estimates of daily potential assimilation to generate a time series of expected LAI that captures the annual phenology for a site. This class provides methods to calculate two time series for LAI:

  • \(L_s\) - the steady state LAI that would be achieved if plants can instantly convert assimilation into leaf area.

  • \(L_r\) - the realised LAI, incorporating lags to capture realistic delays in deploying leaf area.

The class supports different methods for estimating these time series and also provides the Phenology.from_pmodel() method that calculates the required daily assimilation values from a fitted P Model.

Parameters:

  • fapar_limitation – A FaparLimitation instance providing estimates of maximum annual fAPAR and LAI.

  • daily_potential_assimilation – A daily time series of potential assimilation (\(A_0\), mol m-2 day-1).

  • datetimes – The datetimes of the observations of \(A_0\), used to match observations to the appropriate annual values.

  • method – A string selecting the method to be used to estimate LAI.

classmethod from_pmodel(pmodel: PModelABC, fapar_limitation: FaparLimitation, datetimes: ndarray[tuple[Any, ...], dtype[datetime64]] | None = None, gpp_penalty_factor: ndarray[tuple[Any, ...], dtype[floating]] | None = None, **kwargs: Any) → Phenology

Calculate daily phenology from a P Model and other inputs.

TBD.

Parameters:

  • pmodel – A pyrealm.pmodel.pmodel.PModel or pyrealm.pmodel.pmodel.SubdailyPModel instance, fitted with fapar fixed at one.

  • fapar_limitation – A FaparLimitation object providing the maximum annual LAI and fAPAR.

  • datetimes – An array giving the datetimes of observations.

  • gpp_penalty_factor – A GPP penalty factor.

  • **kwargs – Additional arguments.

fapar_limitation

The annual maximum fAPAR and LAI data used in the model.

class PhenologyMethodABC(fapar_limitation: FaparLimitation, daily_potential_assimilation: ndarray[tuple[Any, ...], dtype[floating]], datetimes: ndarray[tuple[Any, ...], dtype[datetime64]], year_index: ndarray[tuple[Any, ...], dtype[int64]], **kwargs: Any)

Abstract base class for implementations of LAI phenology calculation.

Experimental

Be aware that PhenologyMethodABC is an experimental feature and the API and any calculated values may change between major releases.

attrs: tuple[tuple[str, str], ...]

A tuple of attributes to be reported in the Phenology.summarize() output for the method.

method: str

A short method name used to identify the class in PHENOLOGY_METHOD_CLASS_REGISTRY.

requires: tuple[str, ...]

A tuple of any additional variables required to create a Phenology instance using the method.

class PhenologyMethodZhou(fapar_limitation: FaparLimitation, daily_potential_assimilation: ndarray[tuple[Any, ...], dtype[floating]], datetimes: ndarray[tuple[Any, ...], dtype[datetime64]], year_index: ndarray[tuple[Any, ...], dtype[int64]], **kwargs: Any)

Calculation of phenology following (Zhou et al., 2025).

Experimental

Be aware that PhenologyMethodZhou is an experimental feature and the API and any calculated values may change between major releases.

This phenology method class implements the calculation of daily leaf area index from daily assimilation (\(A_0\), mol m-2 day) following the method of (Zhou et al., 2025). The method requires annual values of fAPAR limitation calculated following the method of (Cai et al., 2025) (see FaparLimitationMethodCai).

  • The fractional allocation of GPP to leaf area index (\(m\)) is calculated from the annual maximum LAI and fAPAR as:

    \[m = \frac{\sigma G L_{max}}{A_{0} f_{APAR_{max}}}\]

    Here, $G$ and $A_0$ are the annual estimates of growing season length and annual total potential GPP used to calculate the annual maximum fAPAR and LAI. The coefficient \(\sigma\) is a factor that accounts for the speed with which seasonal canopy growth and senescence patterns affect allocation to LAI (see cai_sigma).

  • The daily allocation of GPP to leaf area index (\(\mu\)) is calculated simply as \(\mu = m A_{d0}\), where \(A_{d0}\) is total daily potential assimilation.

  • The steady state daily estimate of LAI is then calculated as:

    \[L_s = \min \left\{ \mu + \left(\frac{1}{k}\right) W_0 \left[ -k \mu \exp (-k \mu)\right], \text{LAI}_{\text{max}} \right\},\]

    where \(k\) is the light extinction coefficient, \(W_0\) is the principal branch of the Lambert W function. The value is constrained at an upper bound using the estimated annual maximum LAI from the calculate of fAPAR limitation.

  • The realised daily LAI is estimated using an exponential weighted average with a weight \(alpha\) reflecting the time taken for plants to produce leaf area. The value of alpha is set in the phenology constants used to calculate the initial annual fAPAR limitation (zhou_alpha)

attrs: tuple[tuple[str, str], ...] = (('steady_state_lai', '-'), ('realised_lai', '-'))

A tuple of attributes to be reported in the Phenology.summarize() output for the method.

lai_to_gpp_ratio_m: ndarray[tuple[Any, ...], dtype[floating]]

The steady state ratio of leaf area index to potential GPP (\(m\))

method: str = 'zhou'

A short method name used to identify the class in PHENOLOGY_METHOD_CLASS_REGISTRY.

requires: tuple[str, ...] = ()

A tuple of any additional variables required to create a Phenology instance using the method.

class PhenologyMethodZhu(fapar_limitation: FaparLimitation, daily_potential_assimilation: ndarray[tuple[Any, ...], dtype[floating]], datetimes: ndarray[tuple[Any, ...], dtype[datetime64]], year_index: ndarray[tuple[Any, ...], dtype[int64]], aet_pet_ratio: ndarray[tuple[Any, ...], dtype[floating]], spinup_length: int, **kwargs: Any)

Calculation of phenology following (Zhu et al., 2026).

Experimental

Be aware that PhenologyMethodZhu is an experimental feature and the API and any calculated values may change between major releases.

This phenology method class implements the calculation of daily leaf area index from daily assimilation (\(A_0\), mol m-2 day) following the method of (Zhu et al., 2026) The method requires annual values of fAPAR limitation calculated following the method of (Zhu et al., 2026) (see FaparLimitationMethodZhu).

  • The fractional allocation of GPP to LAI (\(m\)) is calculated as the estimated maximum LAI for a year divided by a quantile value (default 0.95) from the distribution of daily assimilation values for that year. The quantile is set when calculating maximum fAPAR via the PhenologyConstNew.zhu_A0_quantile attribute.

  • The steady state daily estimate of LAI is then calculated simply as:

    \[L_s = m A_0.\]

  • The realised daily LAI is estimated using a weighted average over the preceding N days of estimates of steady state LAI. The value of N varies with the climatological estimate of the AET/PET ratio for a site, which must be provided to use this method.

    For the first N observations, this weighted average is calculated over fewer observations - for the first observation it is simply the observed value. In order to avoid this, the method includes an optional spin up that can be used to condition the first values of preceding observation. The conditioning data is taken from the first year of the time series under the assumption that this will be a good approximation for the preceding conditions.

attrs: tuple[tuple[str, str], ...] = (('steady_state_lai', '-'), ('realised_lai', '-'))

A tuple of attributes to be reported in the Phenology.summarize() output for the method.

lai_to_gpp_ratio_m: ndarray[tuple[Any, ...], dtype[floating]]

Annual values for the fractional allocation of potential GPP to leaf area index (\(m\)).

method: str = 'zhu'

A short method name used to identify the class in PHENOLOGY_METHOD_CLASS_REGISTRY.

requires: tuple[str, ...] = ('aet_pet_ratio',)

A tuple of any additional variables required to create a Phenology instance using the method.

PHENOLOGY_METHOD_CLASS_REGISTRY: dict[str, type[PhenologyMethodABC]] = {'zhou': <class 'pyrealm.phenology.phenology.PhenologyMethodZhou'>, 'zhu': <class 'pyrealm.phenology.phenology.PhenologyMethodZhu'>}

A registry for classes implementing different methods for estimating phenology.

Different implementations of the calculation of LAI phenology must all be subclasses of the PhenologyMethodABC abstract base class. This dictionary is used as a registry for defined subclasses and a method name is used to retrieve a particular implementation from this registry. For example:

zhou_phenology_class = PHENOLOGY_METHOD_CLASS_REGISTRY['zhou']