AI-driven prediction of plant physiological traits in hemp using UAV-based multispectral imagery

1 Introduction

Industrial hemp (Cannabis sativa L.) is a versatile bioeconomy crop with whole-plant valorization, including bast fiber and woody hurd for textiles, paper, biocomposites, and construction materials; seed for protein-rich meal and edible oil; and lignocellulosic biomass for bioenergy and bioproducts (De C Lambert et al., 2024; Moscariello et al., 2021). Agronomically, rapid early growth, deep rooting, and dense canopy closure can improve soil structure, suppress weeds, and confer rotational benefits, positioning hemp as an ecologically beneficial component of sustainable systems (Piotrowski and Carus, 2011; Adesina et al., 2020). The crop has been highlighted for strong CO₂ uptake and long-term carbon sequestration potential with modeling estimates suggesting soil organic carbon gains on the order of ∼25.8 t ha^-1 over a 100-year horizon under appropriate management (Muttil et al., 2024; Shen et al., 2022) and it has also been investigated for soil remediation and broader CO₂ bio-sequestration services (Adesina et al., 2020). These environmental attributes align with multiple United Nations Sustainable Development Goals and rural bioeconomy development (Malabadi et al., 2024). Realizing this potential at scale requires optimized nitrogen (N) management to sustain yield and quality while minimizing environmental losses (Kaur et al., 2025). This need highlights the urgency of developing rapid and spatially explicit N diagnostics that support both productivity and environmental sustainability (Tian et al., 2024).

Efficient N management is crucial for sustainable crop production and environmental protection, particularly in emerging industrial crops such as C. sativa (hemp), which is increasingly valued for its fiber, oil, and biomass applications (Antoszewski et al., 2024). In hemp, N availability strongly influences growth, yield, and physiological performance. For example, increasing N fertilization from 50 to 100 kg N ha^-1 increased stem yield by approximately 28%, with smaller gains in inflorescence (≈17%) and seed yield (≈4%) (Campiglia et al., 2017). Similarly, biomass yield has been reported to increase from 1.67 to 4.21 Mg ha^-1 (≈152% increase) when N supply increased from 0 to 200 kg N ha^-1 under temperate conditions (Aubin et al., 2015). In field experiments, seed yield has reached 1,640–2,510 kg ha^-1 at 120–240 kg N ha^-1 (Podder et al., 2024). Nitrogen also affects plant growth and morphology, with plant height increasing from 1.66 to 1.76 m (≈6% increase) and inflorescence length from 66.2 to 82.9 cm (≈25% increase) under N rates of 240 kg ha^-1 compared to unfertilized treatment (Papastylianou et al., 2018). Conventional N diagnostics, such as tissue sampling and laboratory analysis, are destructive, labor-intensive, and lack the spatial resolution required for precision agriculture and variable-rate decision-making (Cheng et al., 2025). Accordingly, high-throughput, non-destructive sensing approaches are a priority for both research and on-farm decision support, and recent reviews underscore the need for standardized, scalable workflows that operate at the field scale and across phenological stages (Yang et al., 2025b).

Recent advances in uncrewed aerial vehicle (UAV) multispectral imaging enable rapid, non-invasive monitoring of canopy traits across field scales with high resolution and revisit frequency (Glenn et al., 2008; Zhou et al., 2023; Zhang et al., 2025a). Vegetation indices (VIs) derived from UAV imagery, including Normalized Difference Vegetation Index (NDVI) and Normalized Difference Red Edge Index (NDRE), have demonstrated strong potential for predicting N status and biomass in crops such as maize and wheat (Cheng et al., 2025; Torres-Sánchez et al., 2014; Jorge et al., 2019; Xie et al., 2023). These results motivate analogous field-scale, UAV-based N studies in hemp production systems.

Recent studies in hemp have primarily focused on structural and compositional traits. For example, UAV multispectral imaging combined with machine learning has been used to estimate plant height, N uptake, and above-ground biomass (Singh et al., 2025), while PROSAIL inversion approaches have retrieved canopy traits such as leaf area index (LAI) and chlorophyll content (Impollonia et al., 2022). Hyperspectral imaging has also been applied for cultivar and growth stage classification (Lu et al., 2022). However, these studies largely emphasize structural or biochemical traits, with limited attention to canopy-scale physiological processes. Consequently, the quantitative prediction of physiological traits such as photosynthesis, stomatal conductance, and transpiration from UAV data remains underexplored in hemp systems. This study addresses this gap by using UAV multispectral data and trait-specific sets of VIs to model SPAD, net photosynthetic rate, stomatal conductance, and transpiration in hemp, thereby extending hemp UAV phenotyping from structural and compositional traits to key physiological functions directly linked to N use and water use.

Physiological traits such as SPAD chlorophyll meter readings (SPAD), net photosynthesis rate (Pn), stomatal conductance (Gs), and transpiration rate (E) are directly linked to plant N metabolism and yield formation (Jones, 2014; Flexas et al., 2012; Xu et al., 2021). Leaf chlorophyll content is closely related to leaf N and photosynthetic capacity, making SPAD-type measures informative for N modeling across phenological stages and genotypes. Classic and recent syntheses demonstrate robust relationships between SPAD, chlorophyll, and photosynthesis that can provide a physiological basis for interpreting spectral signals and improving the biological relevance of UAV-based models (Zhang et al., 2022). Studies also show that integrating Vegetation Indices (VIs) with physiological variables can enhance N prediction accuracy compared to spectral data alone, as demonstrated by an increase in model performance (R ²) from 0.813 to 0.874 and a corresponding reduction in prediction error (Zhang et al., 2025b). Previous studies have demonstrated strong linkages between vegetation indices and plant physiological traits, where indices such as NDVI, NDRE, GNDVI, CIgreen, and CIRE have been successfully used to estimate stomatal conductance, photosynthetic rate, and chlorophyll content in crops including switchgrass (Xu et al., 2021), rice (Wu et al., 2023a), citrus (Xie et al., 2023), maize and soybean (Zhang, 2023), with reported model accuracies ranging from moderate (R ² ≈ 0.29) to high (R ² > 0.90), highlighting the capability of spectral data to capture underlying physiological processes (Xu et al., 2021; Wu et al., 2023b; Xie et al., 2023; Zhang, 2023). Nevertheless, integrative UAV–physiology approaches remain rare in hemp, particularly under diverse N management regimes and across seasons and cultivars.

Machine learning algorithms, including random forest (RF), support vector regression (SVR), partial least squares regression (PLSR), and deep neural networks (DNN), have now been widely applied to estimate leaf N content and related N-status traits across major crops (Zhang et al., 2025b; Hu et al., 2025). In maize, UAV multispectral imagery combined with RF and SVR has predicted canopy N weight with validation R ² values up to 0.85 and RMSE of 4.5 g N m^-2, clearly outperforming linear regression baselines (Lee et al., 2020). In rice, recent deep-learning frameworks that fuse vegetation indices, hyperspectral features, and texture information from UAV imagery have achieved R ² values around 0.87 with RMSE ≈2.6 mg g^-1 for leaf N content (Zhang et al., 2025a), while multi-feature DNN models have improved R ² by ∼40% over conventional machine-learning approaches at key growth stages (Xu et al., 2025). In rice, UAV-based machine learning approaches have also improved N diagnosis, where RF-based models for N nutrition index achieved validation R ² values ranging from 0.61 to 0.79 depending on growth stage (Zha et al., 2020), and similar ML frameworks using UAV imagery have been successfully applied for N estimation (Qiu et al., 2021). In citrus, RF models using UAV-derived spectral features have demonstrated strong predictive performance for canopy N content, achieving R ² values up to 0.90 and outperforming both SVM and linear regression approaches (Prado Osco et al., 2019). Likewise, in winter wheat, machine-learning models integrating UAV multispectral imagery have shown reliable performance for diagnosing N status at field scale (Jiang et al., 2022). A recent meta-analysis synthesizing UAV-based studies across maize, wheat, barley, rice, and oilseed crops further showed that multivariate non-linear methods (e.g., RF, SVM, extreme learning machines) systematically provide stronger effect sizes for N and biomass traits than univariate VI regressions, particularly when models are evaluated on independent validation datasets (Zhang et al., 2024). Comparative studies and meta-analyses across crops generally find that non-linear learners and feature-selection pipelines outperform simple univariate VI regressions, especially when trained on multi-source features and evaluated using independent validation (Zhang et al., 2024; Su et al., 2024; Wang et al., 2024). Despite these advances, applications of UAV-based modeling in hemp remain limited. For instance, Singh et al. (2025) used UAV multispectral imagery combined with machine learning to estimate plant height, biomass, and N uptake in hemp. Similarly, Impollonia et al. (2022) applied PROSAIL model inversion with UAV data to retrieve canopy traits such as LAI and chlorophyll content across cultivars and seasons. In another approach, Lu et al. (2022) employed hyperspectral imaging with machine learning for the classification of hemp cultivars, growth stages, and plant organs. However, these studies have largely focused on structural, biochemical, or classification outcomes, with limited emphasis on modeling key physiological traits such as Pn, Gs, and E specifically using benchmarked low-code analytics platforms for systematic model selection.

Artificial intelligence is increasingly used to turn agricultural sensors and field data into timely predictions while reducing the need for heavy Python programming, which can be a practical barrier for many researchers (Bayar et al., 2025). In agriculture, AutoML (Automated Machine Learning) and low or no code tools automate major steps such as data preparation, feature preprocessing, feature selection, model selection, and hyperparameter optimization, which speeds analysis and reduces manual workload (Shoaib et al., 2023; Sarangpure et al., 2023). Reviews and applications show that these AI workflows help democratize machine learning by making advanced analytics more accessible to non-programmers (Tamayo-Vera et al., 2024; Elbasi et al., 2024). In this study, SAS Viya is used as a low-code environment with an integrated visual interface for building end-to-end machine learning pipelines, including autotuning and automated model comparison (SAS, 2024a; SAS, 2023; SAS, 2024b). Specifically, SAS Viya provides a graphical, workflow-based interface that allows users to develop and compare models using a point-and-click approach, reducing the need for extensive programming while maintaining access to a wide range of standard statistical and machine learning algorithms (e.g., linear regression, generalized linear models, random forests, gradient boosting, and neural networks). This facilitates rapid experimentation, consistent model evaluation, and reproducible workflows, particularly in multidisciplinary research settings where coding expertise may be limited. Recent work on no- or low-code AutoML similarly emphasizes that these platforms enable users with limited coding experience to develop competitive models by automating feature engineering, preprocessing, and tuning (Saranya et al., 2025). However, these platforms do not inherently guarantee improved predictive performance over traditional modeling approaches, and their effectiveness depends on data characteristics and model–trait relationships.

To address the limited integration of UAV-derived vegetation indices with physiological traits and the lack of robust, generalizable approaches for N status prediction in hemp, this study integrates UAV multispectral VIs with ground-measured physiological data across two growing seasons. Specifically, we (i) use UAV multispectral data to calculate trait-specific VI sets to model SPAD, Pn, Gs, and E, and (ii) evaluate the ability of UAV-derived VIs to predict each of these physiological traits independently under contrasting N management strategies. Using SAS Viya, multiple models (e.g., RF, neural networks, regression) are compared for performance and validated through cross-year and cross-cultivar testing across fiber-type and dual-purpose hemp varieties. The goals of this study are to (i) evaluate the ability of UAV-derived vegetation indices to estimate key physiological traits, (ii) identify robust modeling methods for hemp N prediction, and (iii) assess model generalizability across years, varieties, and fertilizer approaches by training models on one subset of data (year or cultivar) and evaluate their predictive performance on independent datasets from different years and hemp types, using metrics such as R ², RMSE, and MAE. By integrating high-throughput UAV-based sensing, physiological ground truth, and AI model selection, this work aims to develop a scalable precision N-monitoring framework to support sustainable hemp production in both research and commercial contexts.

2 Materials and methods

2.1 Experimental design

This study was conducted at the University Farm of North Carolina A&T State University, Greensboro, NC (36°07′N, 79°73′W), during the 2024 and 2025 crop seasons on hemp (C. sativa L.). Dual-purpose (grown for seed and fiber) hemp (C. sativa) varieties the ‘Jurassic’ and ‘Canda’ were evaluated in a randomized complete block design with 4 N application rates (0, 80, 160, and 240 kg N acre⁻¹) and three replications as blocks during 2024 (Figure 1A). The sowing was conducted on 2 May 2024, and emergence occurred 5–6 days after sowing (May 7–8, 2024). Before planting, a base dose of 80 kg N acre⁻¹ was applied to the 80, 160, and 240 kg N treatments, and the remaining N for the 160 and 240 kg N treatments was applied later in the season as a top dressing. Each of the 24 plots contained four rows with 2 ft row spacing, and the plot size was 5 ft × 10 ft (50 ft²). Weed management consisted of weekly manual weeding throughout the growing season. Several light irrigation applications were done using the rain reel system to support early vegetation growth.

Fiber-type hemp varieties, i.e., ‘Yuma’ and ‘Sama’ were evaluated in a randomized complete block design with a split-plot arrangement, with varieties as main plots and 3 N fertility treatments as subplots [control (0 N), inorganic urea, integrated urea plus Nature Safe], and four replications in 2025 (Figure 1B). Sowing was carried out on 24 April 2025, and germination occurred 5–6 days after sowing (April 29–30, 2025). Fertilized treatments targeted 180 kg N ha^-1, with urea applied as approximately 75 percent pre-plant and 25 percent at 1 month after planting; in the integrated treatment, NatureSafe supplied the organic portion at pre-plant, and urea supplied the remainder as a split. Each subplot measured 7 ft × 18 ft (126 ft²; 11.70 m²). Irrigation in 2025 consisted of a single application on 29 April 2025 (approximately 0.10 inch), applied using a rain-reel system to ensure uniform emergence. Weed management consisted of weekly manual weeding throughout the season.

2.2 UAV and multispectral data collection

Multispectral image data were acquired using a WingtraOne GEN II uncrewed aerial vehicle (Wingtra AG, Zurich, Switzerland) equipped with a MicaSense RedEdge-P multispectral camera (MicaSense, Seattle, WA, USA) (Wingtra AG, 2025; MicaSense, 2025c) (Figure 2). The camera was radiometrically calibrated before and after each mission using a calibrated reflectance panel in accordance with the manufacturer’s guidance (MicaSense, 2025a). Flights were planned around local noon to minimize illumination variability, with 80% forward and side overlap (WingtraPilot, 2025). Data collection was scheduled weekly; five flights were completed in 2024, starting 07/16/2024, and eight flights were completed in 2025, starting 06/04/2025, to capture crop growth dynamics across N treatments. All missions were flown at 60 m above ground level to maintain a consistent ground sampling distance. The RedEdge-P recorded six bands comprising blue, green, red, red-edge, near-infrared, and panchromatic (MicaSense, 2025b). Flight plans were preprogrammed in WingtraPilot to ensure complete coverage of the experimental area and sufficient image overlap for accurate mosaicking and processing [WingtraPilot user guide].

2.3 Image processing

Raw multispectral images were first imported into WingtraHub for geotag refinement. Post-processed kinematic geotagging was performed in WingtraHub by pairing the onboard rover observations with RINEX files from a NOAA Continuously Operating Reference Station located about 30 miles from the farm location. The corrected geotags were written to the image EXIF and exported as a calibrated image set for mapping (WingtraHub, 2025; NOAA CORS, 2025). Photogrammetry and reflectance map creation were performed in Pix4D using the standard multispectral workflow. In simple terms, Pix4D detects common features in overlapping photos, aligns them, and stitches them into a single map. Band alignment for the multispectral camera was enabled, and then per-band reflectance orthomosaics were generated at the native ground resolution for each flight date. For plot-level analysis, the orthomosaics were clipped to the experimental boundary and intersected with plot polygons. We extracted zonal statistics (primarily plot-mean values) for each vegetation index and each date. Basic quality control included visual screening to remove any blurred or poorly exposed frames before final processing.

Vegetation indices were calculated in Pix4D’s index calculator so that band math and radiometric handling remained consistent across all flight dates (Figure 3). We computed the 33 vegetation indices listed in Table 1, including broadband greenness indices (e.g., NDVI, OSAVI, SAVI), green- and blue-based chlorophyll proxies (GNDVI, BNDVI, SIPI2), red-edge chlorophyll indices (NDRE, LCI, CIG, Clred-edge, CLRE), structure- and vigor-oriented indices (MGRVI, RGBVI, VARI, TTVI, WDRVI, WDVI, NGRDI, NRVI), chlorophyll-absorption indices (MCARI, MCARI2, TCARI, TGI), and additional spectral metrics such as ARVI, EVI, NDWI, ExG, GLI, CVI, PSRI, RDVI, CCCI, and ICI. NDVI was included as a general proxy for canopy greenness and vigor, with the caveat that it can saturate at high biomass (Tucker, 1979; Huete et al., 2002). To improve sensitivity under dense canopies and to link more directly with chlorophyll and leaf N, we used red-edge indices such as Normalized Difference Red Edge (NDRE) and the Leaf Chlorophyll Index (LCI) alongside related chlorophyll indices (Clevers and Gitelson, 2013; Datt, 1999). Green Normalized Difference Vegetation Index (GNDVI), which leverages the green band, is often more sensitive to chlorophyll across growth stages (Gitelson et al., 1996), while Enhanced Vegetation Index (EVI) was used to reduce soil and background effects and extend dynamic range (Huete et al., 2002). For visible-band options, Visible Atmospherically Resistant Index (VARI) and Triangular Greenness Index (TGI) provide robust responses to canopy greenness and chlorophyll under changing illumination conditions (Gitelson et al., 2002; Hunt et al., 2013). To complement these, Blue Normalized Vegetation Index (BNDVI) offers a blue-based variant that can be informative at early stages or when red reflectance is high (Gordillo-Salinas et al., 2021), and Modified Chlorophyll Absorption in Reflectance Index (MCARI) and related chlorophyll-absorption indices were included to limit background and structural effects when estimating pigment content (Daughtry et al., 2000). All indices in Table 1 were computed from the reflectance mosaics using Pix4D’s default definitions to maintain consistency and reproducibility across dates (Pix4D SA, 2025). These indices were used as predictors to estimate the physiological traits, i.e., SPAD, net photosynthetic rate (Pn), stomatal conductance (Gs), and transpiration (E).

TABLE 1. This table summarizes various vegetation indices and their corresponding mathematical formulas.

To reduce redundancy and focus on the most informative predictors for each physiological trait, we then performed a simple feature-screening step using the combined dataset. For SPAD, Pn, Gs, and E, we calculated the Spearman’s correlation coefficient between each of these traits and each VI (Figures 4A–D). For each trait, the seven indices with the strongest Spearman’s correlations were retained as predictors in the subsequent machine-learning models. The number of selected indices (7) was determined as a practical, data-driven threshold based on diminishing gains in correlation strength and model performance observed during exploratory analysis, while avoiding unnecessary model complexity. This approach is supported by previous studies demonstrating strong relationships between vegetation indices such as NDVI, NDRE, GNDVI, and chlorophyll-based indices with physiological traits including chlorophyll content, photosynthetic activity, and stomatal conductance (Gitelson et al., 2003b; Sims and Gamon, 2002; Wu et al., 2008). The correlation heatmaps presented in this study are based on Spearman’s ρ, providing a robust visualization of VI–physiology relationships.

2.4 Ground truth data collection

On each imaging date, physiological measurements were collected on the same day to align with the UAV acquisitions. Measurements were taken during a consistent mid-day window to reduce diurnal variability and to minimize dew or water-stress effects (for example, 10:30 to 14:30 local time). The sampling order was rotated among N treatments and blocks to avoid systematic time-of-day bias across treatments.

In each plot, one representative plant was selected from the interior rows. Five sunlit, fully expanded leaves were measured per plant using a SPAD-502Plus meter, avoiding the midrib and apparent damage. The five readings were averaged to obtain a single plot-level SPAD value. SPAD measurements were taken on all experimental plots at each sampling date, so that the number of biological replicates per treatment and date equaled the number of replicate plots in the field design. Measurements were collected on clear or lightly overcast days between approximately [09:00 and 14:00] local time to minimize diurnal variation due to low light or midday depression. In 2024, SPAD data were not available for the first two scheduled measurement dates because of instrument constraints. In addition, the cultivar ‘Canada’ was harvested shortly after the second 2024 measurement, and therefore no SPAD readings exist for that cultivar at later 2024 dates. All planned SPAD measurements were successfully collected in 2025, providing a complete time series across all treatments and cultivars for that season.

Plant physiological parameters, i.e., Pn, Gs, and E were measured with a LI-6800 Portable Photosynthesis System (LI-COR Biosciences, Lincoln, NE, USA). To align with the SPAD sampling, one representative plant per plot was selected from the interior rows, and measurements were taken on the topmost fully expanded, sunlit leaf. The chamber was clamped on the middle portion of the lamina, avoiding the midrib and visible defects. For each plot and sampling date, one stabilized reading was recorded after allowing the leaf to equilibrate to chamber conditions (typically within 4–6 min) One stabilized reading per plot was recorded for each variable (Pn, Gs, and E) at each sampling date, and this value was considered the plot-level observation. In 2024, measurements were collected from all plots at each sampling date, resulting in 24 observations on 16 July 2024, and 31 July 2024, and 12 observations on August 12, August 22, and 27 August 2024. In 2025, measurements were collected from all plots across eight sampling dates, with 24 observations recorded at each date. These measurements were carried out on days with stable weather conditions, usually between [09:30 and 14:30] under non-limiting light, to minimize confounding effects of rapidly changing irradiance or leaf temperature. The LI-6800 reference conditions (CO₂ concentration, flow rate, leaf temperature, and relative humidity) were set to approximate ambient conditions at the time of measurement and were kept constant during each measurement sequence within a day.

2.5 Data analysis and modeling

For each flight date, plot-level mean vegetation index values were exported from Pix4D and joined, by plot and date, to the ground-truth table (Figure 3). Modeling was conducted considering SPAD and the physiological variables, i.e., Pn, Gs, and E as the primary target variables. Models were developed in SAS Viya Model Studio (SAS Institute Inc., Cary, NC, USA) (Figure 5). We used a targeted set of supervised learners (i.e., Linear Regression, Generalized Additive Model (GAM), Neural Network, Generalized Linear Model (GLM), Decision Tree, Model Composer (ensemble), Random Forest, Quantile Regression, Gradient Boosting, and Bayesian Additive Regression Trees (BART)) that are well established for tabular agronomic predictors and vegetation-index inputs (Liakos et al., 2018). These models were selected to provide a broad benchmark across interpretable linear models, semi-parametric models, tree-based models, ensemble methods, and neural networks within the SAS Viya framework, reflecting commonly used approaches in remote sensing and crop modeling studies (Panday et al., 2020; Colorado et al., 2020; Kešelj et al., 2026; Latif, 2019; Gaur and Drewry, 2024; Sobejano-Paz et al., 2020; Li et al., 2022). The purpose was comparative screening rather than assuming that all model classes would be equally suitable for the dataset. Decision tree provides an interpretable baseline that captures threshold responses. Random forests improve stability through bootstrap aggregation and handle nonlinearities and correlated indices reliably (Breiman, 2001; Breiman et al., 2017). Gradient boosting adds strong predictive accuracy on structured tabular data by iteratively fitting residuals (Friedman, 2001). A single feed-forward Neural Network was included as the only deep-learning model to provide a representative benchmark for higher-order nonlinear interactions among vegetation indices, while avoiding over-expansion of the model set given the moderate sample size and tabular structure of the dataset, which are generally less favorable for more complex deep-learning architectures. Feed-forward neural networks allow higher-order interactions among indices when the sample size supports it (Goodfellow et al., 2016). The GLM provides an interpretable linear baseline for effect sizes and signs (Nelder and Wedderburn, 1972). GAMs flexibly model smooth nonlinear relationships while retaining partial-dependence-style interpretability (Hastie and Tibshirani, 1986; Wood, 2017). Quantile regression estimates conditional quantiles of the target for uncertainty-aware decision support (Koenker and Bassett, 1978). These families are recommended in agricultural ML surveys for balancing accuracy, interpretability, and robustness across seasons and sites (Liakos et al., 2018; Benos et al., 2021).

Data were split into 60% training, 30% validation, and 10% test, stratified sampling approach rather than simple random partitioning, to preserve the distribution of varieties and N treatments, and as well as sampling dates. This stratification was implemented in SAS Viya Model Studio by selecting the “Stratify” option in the data partitioning step, ensuring that the partitioning process maintained representative distributions of key categorical variables rather than relying on random assignment. (SAS Institute Inc, 2026). This ensured that observations from different sampling periods were proportionally represented across training, validation, and test sets, maintaining the temporal structure of the dataset and ensure that each partition retained representative combinations of cultivar, nitrogen treatment, and time points, and reduce potential bias associated with data partitioning. SPAD data were not available for the first two scheduled measurement dates due to instrument availability constraints, resulting in slight differences in sample size (n) among physiological measurements. In SAS Viya Model Studio (SAS Institute Inc., Cary, NC, USA), candidate models were trained on the training partition and assessed on the validation partition. Model selection was based on test performance rather than training metrics to minimize the risk of overfitting, particularly given the relatively limited dataset size. Outliers were evaluated through visual inspection of data distributions and scatterplots, and no values were removed unless identified as measurement errors, ensuring that biologically meaningful variability was retained.

A feedforward neural network was implemented with one hidden layer consisting of 50 neurons and a hyperbolic tangent (tanh) activation function. Input variables were standardized using midrange scaling, and model training was performed using automatic optimization with hyperparameter tuning enabled. Early stopping based on validation performance was applied using a stagnation criterion of five iterations to reduce the risk of overfitting (SAS Institute Inc., 2026). Tree-based models included decision tree, random forest, and gradient boosting. The decision tree used variance reduction as the splitting criterion, with binary splits, a maximum depth of 20, and a minimum leaf size of 5. The random forest model was configured with 100 trees and an in-bag sampling proportion of 0.6, with other parameters optimized automatically. Gradient boosting models were implemented with 100 trees, a learning rate of 0.1, a subsampling rate of 0.5, and regularization parameters (L1 = 0, L2 = 1), along with early stopping based on a stagnation criterion of five iterations. Tree-based models were further controlled through structural constraints, and all hyperparameters were optimized automatically using validation-based tuning within the SAS Viya framework to reduce overfitting (SAS Institute Inc., 2026).

In addition, statistical models including linear regression, generalized linear models (GLM), generalized additive models (GAM), and quantile regression were implemented for comparison. Linear and quantile regression models used stepwise selection with the Bayesian Information Criterion (BIC). The GLM was configured with a Poisson distribution and log link function. The GAM employed a boosting-based fitting approach with spline effects defined using 3 degrees of freedom and 10 interior knots, with up to 500 iterations and a learning rate of 0.1. All models were evaluated consistently using the same data partitions and performance metrics (SAS Institute Inc., 2026).

2.6 Model evaluation metrics

Prediction performance for each model–trait combination was summarized using the metrics reported by SAS Viya: coefficient of determination (R ²), average squared error (ASE), mean absolute error (MAE), mean squared logarithmic error (MSLE), and root average squared error (RASE).

For a trait with observed values Yi and model predictions Ŷ for i = 1, … , n, and with y ¯ denoting the sample mean of Yi, the statistics were computed as:

R 2 = 1 − ∑ i ˙ = 1 n y i − y ^ i 2 ∑ i = 1 n y i − y ¯ i 2

A S E = 1 n ∑ i ˙ = 1 n y i − y ^ i 2

M A E = 1 n ∑ i = 1 n y i − y ^ i

M S L E = 1 n ∑ i = 1 n ln y i + 1 − ln y ^ i + 1 2

RASE = A S E

R ² expresses the proportion of variance in the trait explained by the model, whereas ASE and MAE describe the average magnitude of squared and absolute residuals, respectively. MSLE emphasize relative rather than absolute errors and stabilizes variance when trait values span several orders of magnitude. All metrics were computed on the independent training, validation, and test set to provide an unbiased assessment of model generalization (Tables 2–5).

TABLE 2. Performance comparison of multiple models for predicting SPAD across training, validation, and testing datasets. Models evaluated include Linear Regression, Generalized Additive Model (GAM), Neural Network, Generalized Linear Model (GLM), Decision Tree, Model Composer (ensemble), Random Forest, Quantile Regression, Gradient Boosting, and Bayesian Additive Regression Trees (BART). Model performance is reported using the coefficient of determination (R ²), average squared error (ASE), mean absolute error (MAE), mean squared logarithmic error (MSLE), and root average squared error (RASE). The dataset was partitioned into training (n = 137), validation (n = 68), and testing (n = 23) subsets. Higher R ² and lower error metrics indicate better model performance.

TABLE 3. Performance comparison of multiple models for predicting Photosynthesis Rate across training, validation, and testing datasets. Models evaluated include Linear Regression, Generalized Additive Model (GAM), Neural Network, Generalized Linear Model (GLM), Decision Tree, Model Composer (ensemble), Random Forest, Quantile Regression, Gradient Boosting, and Bayesian Additive Regression Trees (BART). Model performance is reported using the coefficient of determination (R ²), average squared error (ASE), mean absolute error (MAE), mean squared logarithmic error (MSLE), and root average squared error (RASE). The dataset was partitioned into training (n = 137), validation (n = 68), and testing (n = 23) subsets. Higher R ² and lower error metrics indicate better model performance.

TABLE 4. Performance comparison of multiple models for predicting Stomatal Conductance across training, validation, and testing datasets. Models evaluated include Linear Regression, Generalized Additive Model (GAM), Neural Network, Generalized Linear Model (GLM), Decision Tree, Model Composer (ensemble), Random Forest, Quantile Regression, Gradient Boosting, and Bayesian Additive Regression Trees (BART). Model performance is reported using the coefficient of determination (R ²), average squared error (ASE), mean absolute error (MAE), mean squared logarithmic error (MSLE), and root average squared error (RASE). The dataset was partitioned into training (n = 137), validation (n = 68), and testing (n = 23) subsets. Higher R ² and lower error metrics indicate better model performance.

TABLE 5. Performance comparison of multiple models for predicting Transpiration Rate across training, validation, and testing datasets. Models evaluated include Linear Regression, Generalized Additive Model (GAM), Neural Network, Generalized Linear Model (GLM), Decision Tree, Model Composer (ensemble), Random Forest, Quantile Regression, Gradient Boosting, and Bayesian Additive Regression Trees (BART). Model performance is reported using the coefficient of determination (R ²), average squared error (ASE), mean absolute error (MAE), mean squared logarithmic error (MSLE), and root average squared error (RASE). The dataset was partitioned into training (n = 137), validation (n = 68), and testing (n = 23) subsets. Higher R ² and lower error metrics indicate better model performance.

For statistical evaluation, observed and predicted values of each physiological trait (SPAD, Pn, Gs, and E) were analyzed separately using one-way analysis of variance (ANOVA) to assess the effect of N treatments within each cultivar at each sampling date. When significant treatment effects were detected (p < 0.05), mean separation was performed using Tukey’s honestly significant difference (HSD) test. The results are presented in the graphs, where different letters indicate significant differences among N treatments for each cultivar and date. This analysis was conducted independently for both measured (observed) and model-predicted datasets to evaluate the consistency of treatment responses captured by the models.

3 Results

3.1 Selection of vegetation indices for trait-specific models

Single-VI screening revealed clear differences in the ability of vegetation indices (VIs) to explain variability in physiological traits. Across all traits, indices derived from the green and red-edge spectral regions showed the strongest relationships with SPAD and Pn, whereas Gs and E were better captured by greenness- and canopy-structure-related indices derived from both visible and near-infrared bands (Figures 4A–D).

For SPAD, indices sensitive to chlorophyll content and canopy optical properties consistently ranked among the top-performing predictors (Figure 4A). Spearman’s correlation coefficients for SPAD were very strong and predominantly positive, ranging from ρ = 0.9078 to 0.9375, with the highest correlations observed for LCI (ρ = 0.9375), GNDVI (ρ = 0.9366), and CIG (ρ = 0.9356). In contrast, NDWI showed a strong negative correlation (ρ = −0.9623), indicating an inverse relationship with canopy water-related signals. For Pn, similar patterns were observed (Figure 4B). Strong positive correlations were found for GNDVI (ρ = 0.8931), CIG (ρ = 0.8930), NDVI (ρ = 0.8302), NDRE (ρ = 0.8044), LCI (ρ = 0.8370), and Clred-edge (ρ = 0.7979), while NDWI showed a strong negative correlation (ρ = −0.8934). These results indicate that spectral signals linked to pigment concentration and canopy optical properties are closely aligned with variability in both SPAD and photosynthetic activity. The consistently positive direction of correlations confirms their sensitivity to chlorophyll-driven processes, whereas the negative NDWI relationship reflects contrasting responses to canopy water content. In contrast, Gs showed stronger relationships with indices related to canopy greenness and structural characteristics (Figure 4C). Spearman’s correlations for Gs were moderate to strong, ranging from ρ = 0.7328 to 0.8441, with the highest values observed for MGRVI (ρ = 0.8441), BNDVI (ρ = 0.8316), and NRVI (ρ = 0.8120). Other indices such as VARI (ρ = 0.8118) and WDREVI (ρ = 0.8074) also showed consistent positive relationships, whereas NDVI exhibited a strong negative correlation (ρ = −0.8117). Similar to Gs, indices related to canopy structure and vigor showed the strongest relationships with E (Figure 4D). Spearman’s correlations ranged from moderate to strong (ρ = 0.7519–0.8027), with the highest values observed for TTVI (ρ = 0.8027), WFRAVI (ρ = 0.8026), and NDVI (ρ = 0.8024). Additional indices such as BNDVI (ρ = 0.7839), MGRVI (ρ = 0.7632), and VARI (ρ = 0.7519) also showed positive relationships, whereas NRVI exhibited a strong negative correlation (ρ = −0.8025).

Overall, Spearman’s rank correlation analysis (Figures 4A–D) ensured that only the most informative and physiologically relevant indices were selected for trait-specific model development. The consistently strong positive correlations for SPAD and Pn, compared to moderate relationships for Gs and E, highlight differences in spectral sensitivity among physiological traits and support their subsequent use in predictive modeling frameworks.

3.2 Performance of SPAD prediction models

Using the seven selected SPAD-related VIs within the multi-model framework (including linear, generalized additive, generalized linear, quantile, tree-based, ensemble, and neural-network approaches; Table 2), the SPAD models showed high predictive accuracy on the independent test set (n = 23; Table 2): the best approaches explained more than 88% of the variance (R ² up to 0.88; Table 2) and kept root errors below about 3 SPAD units, meaning predictions typically deviated from measured values by only ∼2 units (Table 2). Across the training (n = 137) and validation (n = 68) datasets, R ² values ranged from 0.68 to 0.96 and 0.61 to 0.78, respectively (Table 2). The lowest training performance was observed for the Neural Network (R ² = 0.68), whereas the highest values were obtained for ensemble approaches such as Model Composer (R ² = 0.96) and Gradient Boosting (R ² = 0.95; Table 2). In the validation dataset, most models showed comparable performance (R ² ≈ 0.73–0.78), including Linear Regression (0.76), Quantile Regression (0.76), GLM (0.74), and Random Forest (0.78; Table 2). The best performer was Linear Regression, with an R ² of 0.883, an ASE of 6.07, MAE of 2.15, MSLE of 0.0062, and corresponding root metric RASE 2.46 (Table 2). Quantile Regression and GLM gave almost identical performance (R ² = 0.880 and 0.878, respectively; Table 2). For Quantile Regression, the test ASE was 6.19, MAE 2.16, MSLE 0.0063, with RASE 2.49; GLM showed ASE 6.30, MAE 2.15, MSLE 0.0066, and RASE 2.51 (Table 2).

Tree-based methods such as Forest and Decision Tree also performed well (R ² = 0.860 and 0.847; Table 2), although with slightly higher prediction errors (ASE = 7.27 and 7.95, respectively; Table 2), while the Neural Network clearly underperformed (R ² = 0.685, ASE 16.32, MAE 2.89; Table 2). Ensemble and machine-learning approaches showed variable performance on the test dataset, with Model Composer (R ² = 0.79), Gradient Boosting (R ² = 0.82), and BART (R ² = 0.82) not exceeding the accuracy of simpler regression models (Table 2). Overall, regression-type models provided the most accurate and parsimonious SPAD predictions across training, validation, and testing datasets (Table 2).

The relationship between measured and Linear Regression–predicted SPAD was strongly linear (Figure 6). Observations were closely aligned along the 1:1 line across the full SPAD range, with minimal deviation at the lower and upper extremes. The 95% confidence band around the regression line was narrow, indicating limited uncertainty in the mean prediction, whereas the wider 95% prediction band encompassed nearly all observations. This visual pattern is consistent with the high test performance of the Linear Regression model reported in Table 2.

3.3 Performance of Pn prediction models

Using the seven selected Pn-specific VI set and modelling framework (Table 3), the models for Pn also showed strong predictive skill on the test set (n = 28; Table 3). The best models explained about 85%–87% of the variance (R ² ≈ 0.85–0.87; Table 3) and had root errors of only ∼2–3 µmol CO₂ m^-2 s^-1, indicating that UAV-derived VIs captured most of the plot-to-plot variation in Pn (Table 3).

Across the training (n = 165) and validation (n = 83) datasets, R ² values ranged from 0.64 to 0.95 and 0.63 to 0.82, respectively (Table 3). The lowest performance was observed for the Neural Network (R ² = 0.64 training; 0.63 validation), whereas the highest values were obtained for Gradient Boosting (R ² = 0.95 training) and Linear Regression (R ² = 0.82 validation; Table 3). Most models showed comparable validation performance (R ² ≈ 0.70–0.82), including GLM (0.74), Random Forest (0.80), and BART (0.80). The highest R ² on the test dataset was obtained with GAM (R ² = 0.87; Table 3), which had ASE 7.39, MAE 2.13, MSLE 0.10, and RASE 2.72. Gradient Boosting performed very similarly (R ² = 0.87; Table 3), with ASE 7.53, MAE 1.96, MSLE 0.12, and RASE 2.74. Linear Regression and Quantile Regression also achieved high test R ² values of 0.86 and 0.85, respectively (Table 3). For Linear Regression, ASE and MAE were 8.07 and 2.28 (MSLE 0.10, RASE 2.84), while Quantile Regression showed ASE 8.24, MAE 2.31, MSLE 0.10, and RASE 2.87 (Table 3). The remaining models had slightly lower performance (R ² ≈ 0.81–0.85; Table 3), while the Neural Network and GLM showed comparatively lower accuracy (R ² ≈ 0.74; ASE ≈14.4; Table 3). Overall, multiple modelling approaches provided robust predictions of Pn, with GAM and Gradient Boosting showing the highest accuracy on the test dataset (Table 3).

The relationship between measured and GAM-predicted photosynthesis rate was also strongly linear (Figure 7). Observations were closely aligned along the 1:1 line across the full range of measured Pn, with only minor deviations at higher values. The 95% confidence band around the regression line was narrow, suggesting limited uncertainty in the mean prediction, whereas the wider 95% prediction band encompassed nearly all observations. No strong systematic bias was observed, confirming the good predictive performance of the GAM model for Pn reported in Table 3.

3.4 Performance of Gs prediction models

Within the same multi-model framework (Table 4), Gs was predicted with consistently good accuracy on the test dataset (n = 28; Table 4), with the best models explaining just over 80% of the variance (R ² ≈ 0.81–0.82; Table 4) and maintaining root errors around 0.11–0.12 mol H₂O m^-2 s^-1 (Table 4).

Across the training (n = 165) and validation (n = 83) datasets, R ² values ranged from 0.72 to 0.92 and 0.68 to 0.81, respectively (Table 4). The lowest performance was observed for the Neural Network (R ² = 0.72 training; 0.68 validation), whereas the highest values were obtained for Gradient Boosting (R ² = 0.92 training) and Model Composer (R ² = 0.81 validation; Table 4). Most models showed comparable validation performance (R ² ≈ 0.70–0.81), including Linear Regression (0.80), Quantile Regression (0.80), and Decision Tree (0.80). The best results on the test dataset were obtained with Quantile Regression (R ² = 0.82; Table 4), which had ASE 0.013, MAE 0.086, MSLE 0.006, and RASE 0.115. Linear Regression performed similarly (R ² = 0.82; Table 4), with ASE 0.013, MAE 0.085, MSLE 0.006, and RASE 0.116. The Decision Tree model also showed strong performance (R ² = 0.81; Table 4), with ASE 0.014 and MAE 0.092, while the Model Composer ensemble was slightly lower (R ² = 0.81; Table 4). Other models showed comparable but slightly lower performance (R ² ≈ 0.75–0.80; Table 4), while the Neural Network and GLM had the lowest predictive accuracy (R ² ≈ 0.70–0.74; Table 4). Overall, multiple modelling approaches provided reliable predictions of Gs, with Quantile Regression and Linear Regression showing the highest accuracy on the test dataset (Table 4).

The relationship between measured and Quantile Regression–predicted stomatal conductance was also highly linear (Figure 8). Observations were closely aligned along the 1:1 line across the full range of measured Gs, with only minor increases in scatter at higher values. The 95% confidence band around the regression line remained narrow, while the broader 95% prediction band encompassed nearly all points. No clear systematic bias was observed, consistent with the strong predictive performance of the Quantile Regression model (R ² ≈ 0.82; Table 4).

3.5 Performance of E prediction models

E was the most difficult trait to predict within the same multi-model framework (Table 5), but several algorithms still provided useful accuracy on the test dataset (n = 28; Table 5). The best models explained about 68%–75% of the variance (R ² ≈ 0.68–0.75; Table 5) and yielded very small absolute and root errors (MAE ≈0.014–0.017; RASE ≈0.017–0.020; Table 5), indicating that predictions captured a substantial portion of variation in E.

Across the training (n = 165) and validation (n = 83) datasets, R ² values ranged from −0.01 to 0.92 and −0.01 to 0.79, respectively (Table 5). The lowest performance was observed for GLM (R ² ≈ −0.01), whereas the highest values were obtained for BART (R ² = 0.92 training) and Model Composer (R ² = 0.79 validation; Table 5). Most models showed moderate validation performance (R ² ≈ 0.69–0.79), including Linear Regression (0.76), Quantile Regression (0.76), and Gradient Boosting (0.77). The best performance on the test dataset was obtained with Quantile Regression (R ² = 0.75; Table 5), which had ASE 0.0003, MAE 0.014, MSLE 0.0003, and RASE 0.017. The Model Composer ensemble (R ² = 0.74; Table 5) and GAM (R ² = 0.74; Table 5) performed similarly, with ASE ≈0.0003 and MAE ≈0.015. Linear Regression also showed comparable performance (R ² = 0.73; Table 5), with ASE 0.0003, MAE 0.015, and RASE 0.018. Other models showed slightly lower performance (R ² ≈ 0.66–0.69; Table 5), including Neural Network and Random Forest, while GLM produced near-zero or negative R ² values (≈−0.02; Table 5), indicating poor predictive performance. Overall, although predictive accuracy for E was lower than for SPAD, Pn, and Gs, several models still provided reliable estimates with low absolute error (Table 5).

The relationship between measured and Quantile Regression–predicted transpiration rate was also strongly linear (Figure 9). Observations were closely aligned along the 1:1 line across the full range of measured E, with slightly greater dispersion at higher values. The 95% confidence band around the regression line was narrow, indicating low uncertainty in the mean prediction, while the wider 95% prediction band encompassed nearly all observations. No clear systematic bias was observed, consistent with the predictive performance of the Quantile Regression model (R ² ≈ 0.75; Table 5).

3.6 Measured and predicted physiological parameters across the N treatments

Physiological traits varied in response to N treatments across cultivars and growing seasons, with strong agreement between observed and predicted values. Chlorophyll content (SPAD) increased consistently with N application in both years (Figures 10a–f). In 2024 (Jurassic cultivar), SPAD values were higher under increased N rates, particularly at early sampling dates, but differences among treatments diminished toward later stages, indicating a saturation effect (Figure 10a for observed; Figure 10d for predicted; 2024 Jurassic). In 2025, both Yuma and Sama cultivars exhibited stronger temporal dynamics, with SPAD increasing from early growth stages, peaking around late June, and declining thereafter (Figures 10b,c for observed; Figures 10e,f for predicted; 2025 Yuma and Sama). Higher N treatments (T2 and T3) consistently resulted in greater SPAD values than T1, with clearer separation among treatments compared to 2024. The Yuma cultivar showed greater variability and stronger responsiveness, whereas Sama displayed a more stable and gradual pattern (Figures 10d–f). Predicted SPAD values closely followed observed trends, with slight smoothing at peak values (Figures 10d–f).

Pn also increased with N application, although responses varied across cultivars and sampling dates (Figures 11a–h). In 2024, higher N rates generally enhanced Pn, particularly at early growth stages, but treatment differences became less consistent later in the season (Figures 11a,c for observed; Figures 11d,f for predicted; 2024 Jurassic and Canda). In contrast, 2025 showed stronger N effects, with Pn increasing toward a peak around late June before declining (Figures 11e,g for observed; Figures 11f,h for predicted; 2025 Yuma and Sama). Higher N treatments consistently produced greater Pn values, especially during peak growth. The Yuma cultivar exhibited higher Pn and greater variability compared to Sama, indicating stronger sensitivity to N and environmental conditions. Predicted Pn values captured these trends well, maintaining treatment ranking with minor underestimation at higher values (Figures 11d,f,h).

In contrast, Gs and E did not show significant differences across N treatments in either observed or predicted datasets (Figures 12a–h; Figures 13a–h). Similar statistical groupings across treatments indicate that these parameters were largely unaffected by N levels. Both Gs and E remained relatively stable across sampling dates, cultivars, and years (Figures 12a,c,e,g; Figures 13a,c,e,g for observed; Figures 12b,d,f,h; Figures 13b,d,f,h for predicted; representing 2024 Jurassic and Canda, and 2025 Yuma and Sama), suggesting that stomatal behavior and water flux were less responsive to N variability compared to photosynthetic capacity and chlorophyll content. Predicted values closely matched observed data, reflecting the limited variability in these traits (Figures 12b,d,f,h; Figures 13b,d,f,h).

4 Discussion

This study shows that UAV-derived multispectral vegetation indices (VIs) can recover key physiological traits of field-grown hemp with high accuracy, with the best models achieving test-set R ² values of 0.88 for SPAD, 0.87 for Pn, 0.82 for Gs, and 0.75 for E. Using trait-specific groups of seven VIs, we explained nearly 90% of the variation in SPAD and Pn, over 80% of the variation in Gs, and about 75% of the variation in E on an independent test set. These performance levels are comparable to, or better than, those reported for other crops using hyperspectral sensors or combined optical–thermal approaches (Furbank et al., 2021; Sobejano-Paz et al., 2020; Parkash et al., 2025). To our best knowledge, this is one of the first studies to showcase that UAV multispectral data alone can be used to predict SPAD, Pn, Gs, and E simultaneously in hemp.

4.1 UAV multispectral phenotyping in hemp

Recent UAV-based studies across multiple crops have demonstrated the capability of multispectral imagery combined with machine-learning approaches to estimate key physiological traits, including chlorophyll content (SPAD), photosynthetic rate (Pn), and stomatal conductance (Gs). For example, studies in wheat, rice, and potato have shown that models such as random forest, gradient boosting, neural networks, and partial least squares regression can reliably capture these physiological responses across growth stages, often with improved performance when combining spectral, structural, and temporal features (Wu et al., 2023a; Wu et al., 2023b; Yin et al., 2023; Zhang et al., 2025a). These findings highlight the growing potential of UAV-based approaches to move beyond structural trait estimation toward functional, physiology-based phenotyping. Impollonia et al. (2022) used PROSAIL inversion of UAV multispectral images to estimate LAI and leaf chlorophyll content in hemp, while Singh et al. (2025) coupled UAV multispectral imagery with machine-learning models to predict plant height, leaf N uptake, and above-ground biomass across N fertilization treatments. Serrano et al. (2025) conducted a 3-year N-rate experiment in floral hemp (‘Wife’) and used UAV-based multispectral imaging to quantify biomass responses to 6 N application rates (0–280 kg N ha^-1). They showed that higher N rates (e.g., 224 kg N ha^-1) significantly increased above-ground biomass relative to lower rates, and that canopy area, NDVI, and GNDVI measured near harvest were strongly and positively correlated with biomass, demonstrating that simple VIs can sensitively track N-driven differences in crop growth and support N management decisions. However, despite these advances, hemp studies have largely focused on structural and compositional traits captured at limited time points. By contrast, our study uses UAV-derived VIs to predict canopy-scale physiological processes SPAD, Pn, Gs, and E across multiple dates, thereby extending hemp UAV applications from static structural and chlorophyll traits to dynamic gas-exchange behavior and filling a gap in the existing literature.

By linking multispectral imagery directly to SPAD, Pn, Gs, and E, our work advances hemp phenotyping into the physiological realm. This brings hemp closer to the ‘physiology-aware’ phenotyping currently being developed in wheat and cotton, where imaging is used to estimate photosynthetic capacity and gas-exchange behavior (Furbank et al., 2021; Parkash et al., 2025). Using only seven VIs per trait is also practically important. Breeding and management workflows often need to analyze many flights and thousands of plots; calculating dozens of indices per image can be both computationally and analytically burdensome (Xie and Yang, 2020; Ayankojo et al., 2023). Our results show that a small, carefully selected set of indices can be enough to capture most of the spectral information relevant to hemp chlorophyll status and gas exchange. This aligns with recent hemp pipelines that emphasize streamlined extraction of a limited number of robust spectral features for downstream genetic or management analyses (Impollonia et al., 2022; Singh et al., 2025).

4.2 Spectral controls on SPAD and gas-exchange traits

The strong relationships between SPAD and red-edge or green-based indices, including CIG, GNDVI, NDWI, Clred-edge, LCI, and NDRE, align with the long-recognized sensitivity of these indices to chlorophyll and N content (Clevers and Gitelson, 2013; Zhang et al., 2021; Lin et al., 2019). CIgreen and CIred-edge were originally developed as narrow-band chlorophyll indices that remain sensitive at high pigment levels, where classic indices like NDVI tend to saturate (Ciganda et al., 2009; Clevers and Gitelson, 2013; Xie et al., 2018). Spearman’s correlation analysis further confirmed these relationships, with strong positive correlations for GNDVI (ρ = 0.9366), CIG (ρ = 0.9356), and LCI (ρ = 0.9375), and a strong negative correlation for NDWI (ρ = −0.9623), indicating its inverse relationship with chlorophyll content. Physiologically, this strong relationship is expected because SPAD provides a proxy for leaf chlorophyll and, indirectly, leaf N per unit area, so between-plot differences in SPAD largely reflect differences in leaf N status that co-vary with canopy photosynthetic capacity (Clevers and Gitelson, 2013; Lin et al., 2019). Moreover, the red-edge region sits on the steep transition between strong red absorption and high near-infrared reflectance, so small changes in chlorophyll concentration cause relatively large changes in red-edge reflectance, helping CIG capture subtle differences in hemp canopy greenness even at high pigment levels (Clevers and Gitelson, 2013; Xie et al., 2018).

The same set of indices also ranked highly for Pn. This is supported by strong Spearman correlations, where GNDVI (ρ = 0.8931), CIG (ρ = 0.8930), and NDVI (ρ = 0.8302) showed strong positive relationships with Pn, while NDWI exhibited a strong negative correlation (ρ = −0.8934). The negative relationship between NDWI and SPAD and Pn reflects the different physiological and optical mechanisms captured by these indices. NDWI is primarily sensitive to canopy water content and internal leaf structure, driven by water absorption features in the near-infrared region (Gao, 1996; Peñuelas et al., 1993), whereas SPAD and Pn are closely linked to chlorophyll concentration and nitrogen-driven photosynthetic capacity (Evans, 1989; Croft et al., 2017). Under higher nitrogen availability, increases in chlorophyll content and photosynthetic activity are not necessarily accompanied by proportional increases in canopy water signals, resulting in an apparent inverse relationship between NDWI and SPAD/Pn. Additionally, nitrogen-induced changes in leaf anatomy, including mesophyll structure and light scattering properties, may further decouple water-related spectral responses from chlorophyll-driven reflectance signals (Gitelson et al., 2003a). Similar relationships between red-edge/green indices, chlorophyll, and light-saturated photosynthesis have been documented in wheat, rice, and cotton (Furbank et al., 2021; Lin et al., 2019; Parkash et al., 2025). Physiologically, this aligns with expectations: in C₃ crops, variation in photosynthetic capacity often correlates with leaf N and chlorophyll, especially under non-severe stress, so indices that measure pigment pools tend to reflect photosynthetic potential.

Vegetation index saturation, particularly for indices such as NDVI under high canopy density, is a known limitation in optical remote sensing. In this study, this effect was partially mitigated by including vegetation indices derived from red-edge and green spectral bands (e.g., CIred-edge, GNDVI, CIG), which are less prone to saturation under moderate-to-high biomass and chlorophyll conditions (Xie et al., 2018). In addition, the use of multi-temporal observations across sampling dates ensured that model training incorporated variability across growth stages rather than relying solely on late-stage canopy conditions. Furthermore, the modeling framework utilized multiple indices simultaneously, allowing the models to capture complementary spectral information and reduce dependence on any single index that may be affected by saturation.

In contrast, Gs and E were best captured by a slightly different spectral set. For Gs, indices such as MGRVI (ρ = 0.8441), BNDVI (ρ = 0.8316), and NRVI (ρ = 0.8120) showed strong positive correlations, while NDVI showed a negative relationship (ρ = −0.8117). Similarly, for E, NDVI (ρ = 0.8024), TTVI (ρ = 0.8027), and WFRAVI (ρ = 0.8026) showed consistent moderate-to-strong correlations. This shift aligns with studies showing that canopy-scale indicators of structural vigor and water status, such as NDVI, WDRVI, and moisture-sensitive indices, are effective proxies for stomatal regulation and transpiration in cereals and citrus (Sobejano-Paz et al., 2020; Xie et al., 2023; Yang et al., 2025a). Stomatal conductance is influenced by leaf water content, turgor, and boundary-layer conditions, all closely related to canopy density and leaf area, and thus to green–NIR reflectance. Yang et al. (2025b) and others have also emphasized that combining structural VIs with meteorological data is especially effective for predicting Gs and E at the canopy level. The relatively limited response of Gs and E to nitrogen treatments observed in this study is consistent with the physiological regulation of these traits. Unlike chlorophyll content and photosynthetic capacity, which are strongly influenced by nitrogen availability, stomatal conductance and transpiration are primarily controlled by environmental drivers such as vapor pressure deficit, temperature, and radiation, as well as plant hydraulic status (Medlyn et al., 2011; Buckley, 2019).

Overall, the trait-specific VI rankings indicate that red-edge-based chlorophyll indices are particularly effective at identifying pigment pools and photosynthetic capacity, while structural and “vigour” indices provide more insight into traits controlled by stomatal regulation and water movement. This aligns with the common understanding that spectral proxies for chlorophyll, structure, and water content each reflect different aspects of plant function (Xie and Yang, 2020; Lin et al., 2019; Yang et al., 2025b). Vegetation indices derived from similar spectral bands are expected to exhibit multicollinearity due to overlapping spectral sensitivity. However, multicollinearity primarily affects the interpretability of individual model coefficients rather than predictive performance (Dormann et al., 2013). Since this study focused on prediction rather than parameter inference, and given the high redundancy among vegetation indices, model evaluation was based on predictive performance rather than coefficient-level interpretation. To further examine redundancy among vegetation indices, LASSO-based feature selection was performed for each physiological trait. For SPAD and Pn, CIG was consistently identified as the dominant predictor, with additional indices providing minimal improvement in model performance. When CIG was excluded, multiple alternative indices were selected and predictive accuracy remained high, confirming strong redundancy among chlorophyll-related indices. In contrast, for Gs and E, fewer predictors were selected and model performance was comparatively lower, reflecting the more complex physiological regulation of gas-exchange traits. These results indicate that multicollinearity is inherent in UAV-derived spectral datasets and that predictive performance is driven by a dominant underlying spectral signal rather than independent contributions of individual indices. Similar observations have been reported in UAV-based crop studies, where multicollinearity among vegetation indices is expected and predictive models remain robust despite correlated inputs (Derraz et al., 2023).

Although this study focused on predictive performance rather than detailed feature attribution, the contribution of explanatory variables can be interpreted based on the consistent ranking of vegetation indices and their known physiological relevance. Red-edge- and green-based indices (e.g., CIG, GNDVI, NDRE) showed the strongest relationships with SPAD and Pn, reflecting their sensitivity to chlorophyll content and photosynthetic capacity. In contrast, indices associated with canopy structure and water status (e.g., NDVI, WDRVI, NDWI) contributed more to the prediction of Gs and E. While formal variable importance or explainable AI techniques (e.g., SHAP values) were not applied in this study, the observed relationships are consistent with established spectral–physiological linkages and provide a qualitative interpretation of feature contributions. Future work should incorporate quantitative feature importance analyses to further refine these relationships.

4.3 Choice of modelling approach

An interesting result is that relatively simple regression-type models often matched or outperformed more complex machine-learning algorithms. Across both training and validation datasets (Tables 2–5), linear regression, GAM, and quantile regression consistently showed stable performance, indicating that simpler models were sufficient to capture the relationships between VIs and physiological traits. For SPAD, a multiple linear regression based on seven VIs achieved a test R ² of approximately 0.88, slightly surpassing random forests and gradient boosting. For Pn, the best-performing model was a generalized additive model (GAM), while for Gs and E, quantile regression achieved the highest accuracy.

There are several plausible reasons. First, the relationships between the selected VIs and the traits were mostly smooth and monotonic, as shown by the strong measured–predicted relationships and narrow confidence bands in the scatterplots. In this context, flexible yet still parametric methods such as linear regression and GAMs can effectively capture the main signal without the overfitting risk associated with high-capacity learners, especially when sample sizes are moderate (Li et al., 2022; Zhou et al., 2021). Second, the VI sets were intentionally pruned to reduce redundancy, so tree-based models had less opportunity to leverage complex interactions compared to studies using hundreds of spectral bands (Furbank et al., 2021; Parkash et al., 2025). In addition, vegetation indices are engineered spectral features that already summarize physiologically meaningful canopy information, and their relationships with traits such as SPAD and Pn are often approximately monotonic. This explains why simpler models, such as linear and quantile regression, performed comparably to more complex approaches, as the predictive signal was largely captured by a limited number of dominant indices (Peng et al., 2017).

However, several recent studies have shown little improvement over simpler models once the most informative bands or indices are identified (Xie et al., 2023; Parkash et al., 2025). For breeding and agronomic applications, the interpretability, computational efficiency, and ease of deployment of linear or GAM-type models are key advantages. Our findings suggest that, for hemp phenotyping with a multispectral camera, these simpler models are a practical default, with more complex algorithms used only when trait–spectral relationships are clearly nonlinear or the training data cover many environments. The inclusion of multiple modeling approaches in this study was intended as an exploratory benchmarking step within the SAS Viya framework to systematically compare linear, semi-parametric, and machine-learning methods under consistent conditions, rather than to suggest that all model types are necessary for similar tabular remote sensing applications.

The ensemble “Model Composer,” which aggregates predictions from multiple base learners (e.g., linear, tree-based, and boosting models) into a single combined model, showed very high training R ², especially for SPAD, but did not outperform the best individual learners on the test set, indicating mild overfitting. This overfitting behavior can be attributed to the relatively limited dataset size and the use of a small set of highly informative but partially correlated vegetation indices, which may lead high-capacity models to fit noise in the training data. Although strategies such as early stopping, hyperparameter tuning, and regularization were implemented within the SAS Viya framework, their effectiveness is inherently constrained by the amount and variability of available data. Under such conditions, simpler models with stronger structural constraints tend to generalize better, which is consistent with the observed performance in this study. This highlights that ensembles are not guaranteed to beat simpler models in moderately sized phenotyping datasets and should be critically evaluated rather than assumed to be superior. It is important to note that the use of an AutoML framework in this study was intended to facilitate systematic model development and comparison rather than to demonstrate superiority over traditional modeling approaches. The models implemented in SAS Viya (e.g., linear regression, generalized additive models, and machine-learning algorithms) are equivalent to those commonly developed in code-based environments; however, they were applied here through a graphical, point-and-click workflow that standardizes preprocessing, partitioning, and evaluation. This approach enables consistent and reproducible model comparison while reducing the need for extensive manual coding. Accordingly, the objective of this study was to identify robust relationships between UAV-derived vegetation indices and physiological traits under a unified modeling framework, rather than to benchmark AutoML against manually tuned implementations, which was beyond the scope of this work.

4.4 Implications for hemp breeding and crop management

Being able to estimate SPAD, Pn, Gs, and E from UAV imagery opens several opportunities for hemp improvement and management. These traits are central to the crop’s carbon–water economy and are crucial for biomass accumulation and, potentially, cannabinoid yield (Chandra et al., 2011). Incorporating our trait predictions into routine flights would allow breeders to screen large hemp populations for photosynthetic efficiency and stomatal behavior, in addition to plant height or biomass. Similar strategies in wheat and other cereals have sped up selection for canopy photosynthesis and radiation-use efficiency (Wang et al., 2020; Volpato et al., 2021; Cheng et al., 2025).

For growers, the fairly high accuracy of Gs and E models indicates that multispectral UAV data could help monitor water status in hemp fields. Reviews of UAV-based crop water monitoring emphasize the importance of combining VIs with thermal data and weather information to map spatial patterns of stomatal closure and transpiration (Yang et al., 2025b; Yang et al., 2025a). Although we did not include thermal imagery here, our results show that optical VIs alone already capture much of the variation in Gs and E across plots. In practice, spectral models like ours could identify areas with limited water availability or damaged root function, supporting variable-rate irrigation or targeted soil management.

Importantly, most of the high-ranking VIs in our study, CIG, NDRE, GNDVI, NDVI, NDWI, WDRVI, and related indices, can be calculated using the standard green, red, red-edge, and NIR bands available on many commercial agricultural UAVs. This aligns our approach with existing hemp imaging workflows that have used similar bands to estimate LAI, chlorophyll, and biomass (Impollonia et al., 2022; Serrano et al., 2025). As a result, the models developed here could be implemented without specialized hardware, making them appealing for both research and on-farm use.

While the present study demonstrates reliable prediction of physiological traits under varying nitrogen treatments, direct translation of these results into precision nitrogen management tools (e.g., diagnostic thresholds, critical nitrogen dilution curves, or variable-rate fertilization strategies) was beyond the scope of this work. In this study, nitrogen was used as a treatment factor to induce variability in plant responses rather than as a direct prediction target. However, the predicted spatial and temporal variability in SPAD and photosynthetic traits provides a practical basis for identifying zones of differential crop performance within fields. These zones can be used to guide targeted field scouting, soil or tissue sampling, and site-specific management decisions, thereby supporting more informed nitrogen management strategies. In addition, temporal monitoring of these physiological traits may enable early detection of nutrient-related stress before visible symptoms appear, improving the timing of management interventions. Future research should integrate UAV-derived physiological traits with plant nitrogen measurements to establish diagnostic thresholds and enable practical nitrogen management applications in hemp.

4.5 Limitations and future directions

Several limitations must be acknowledged. First, although multiple machine learning models were evaluated, the relatively limited dataset size may constrain model generalizability and increase the risk of overfitting, and therefore the predictive performance should be interpreted with caution as a comparative assessment rather than absolute prediction accuracy. In addition, the models were calibrated and tested in a single environment. Relationships between VIs and physiological traits can change depending on genotype, soil background, illumination geometry, and atmospheric conditions (Xie and Yang, 2020; Xie et al., 2018; Yang et al., 2025b). Multi-year, multi-site datasets covering various hemp cultivars are necessary to evaluate model transferability and, if needed, develop domain-adaptation strategies.

Second, we used a multispectral camera with a limited number of broad bands. Hyperspectral systems capture more detailed absorption features and continuum shapes associated with pigments, water, and structural traits, and have been shown to improve predictions of photosynthetic parameters in several species (Furbank et al., 2021; Parkash et al., 2025; Zhou et al., 2021). Future work could compare our compact VI-based models with hyperspectral or combined optical–thermal approaches to quantify the trade-offs between additional accuracy and operational complexity in hemp phenotyping.

Third, we modeled plot-average traits measured at a single time of day. Gas-exchange traits are dynamic and respond to short-term fluctuations in light, temperature, and vapor pressure deficit. Diurnal and seasonal patterns may therefore not be fully represented in our dataset. Incorporating time-series UAV flights and coupling them with mechanistic or explainable machine-learning models of stomatal conductance (Gaur and Drewry, 2024; Li et al., 2022; Mayanja et al., 2026) would help distinguish genotypic differences from short-term environmental noise.

Finally, although prediction errors were small for most traits, residual uncertainty remains, especially for E. For applications like precision irrigation, acceptable error margins might be tighter than for breeding selection. Combining our spectral models with additional data such as canopy temperature, soil electrical conductivity, or fine-scale topography could further decrease uncertainty and improve robustness in operational settings (Yang et al., 2025b; Xie et al., 2023).

4.6 Conclusions

In summary, this study demonstrates that UAV-based multispectral imaging can effectively predict key physiological traits of hemp under varying nitrogen conditions, with strong predictive performance for SPAD and net photosynthetic rate, and moderate accuracy for stomatal conductance (Gs) and transpiration (E). Trait-specific vegetation indices showed consistent relationships with physiological measurements, where red-edge and green-based indices (e.g., CIG, GNDVI, NDRE) were most effective for chlorophyll-related traits, while structure- and vigor-related indices provided relatively better insights into Gs and E.

Across modeling approaches, simpler models such as linear regression, generalized additive models (GAM), and quantile regression performed comparably to, and in some cases better than, more complex machine-learning methods. This suggests that the relationships between selected VIs and physiological traits were largely smooth and could be captured without high model complexity. However, given the relatively limited dataset size, these results should be interpreted as demonstrating relative model performance rather than definitive predictive capability across broader environments.

Overall, these findings highlight the potential of UAV-based multispectral phenotyping as a practical tool for estimating physiological traits in hemp. The proposed framework provides a basis for integrating UAV-derived indicators into precision nitrogen management and breeding programs, although further validation across multi-site and multi-season datasets is required to confirm model robustness and generalizability.

Generative AI statement

The author(s) declared that generative AI was not used in the creation of this manuscript.

Any alternative text (alt text) provided alongside figures in this article has been generated by Frontiers with the support of artificial intelligence and reasonable efforts have been made to ensure accuracy, including review by the authors wherever possible. If you identify any issues, please contact us.

Publisher’s note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Footnotes

Edited by: Yu Fenghua, Shenyang Agricultural University, China

Reviewed by: Suchitra Mohan Patil, K. J. Somaiya College of Engineering, India

Shengzhao Pei, Northwest A&F University, China