The experiments were carried out during the production cycles of autumn-winter 2017/2018 and 2018/2019 at the La Laguna Experimental Station in the Instituto Nacional de Investigaciones Forestales, Agricolas y Pecuarias (INIFAP), located in Matamoros, Coahuila, Mexico (25° 32’ N, 103° 14’ O and 1150 m above sea level). The soil at the experimental site was a loamy-clayed texture, with a depth greater than 1.8 m. The availability of water was 150 mm m^-1 [13], the organic carbon content was 0.75%, and the pH was 8.14 [14].

The research consisted of the comparison between the behavior of four safflower cultivars harvested in four phenological stages under a randomized complete blocks design with four repetitions in a 4 × 4 factorial arrangement (64 experimental plots). Factor A was the cultivars, and Factor B was the phenological stage. The cultivars were: CD868, Selkino, Gila, and Guayalejo, with the first two being spineless and the other two spiny; all of them were developed by INIFAP. The phenological stages during harvest were: (E): beginning of capitulum formation (E50), capitulum clearly separated from the younger leaves (E55), distinguishable medium and intermediate external bracts (E59), and beginning of flowering (E61) [15]. The harvest was made 76, 84, 92, and 100 days after sowing (das) in E50, E55, E59, and E61, respectively. All cultivars were harvested simultaneously at each harvest date.

Seedbed preparation was performed through disk plough at a depth of 0.30 m, followed by double disking and zero-slope levelling. Nitrogen and phosphorus fertilizer dose was calculated considering the safflower extraction capacity: 250 kg N ha^-1 and 80 kg P₂O₅ ha^-1. The nitrogen source was urea (46% N), and monoammonium phosphate (52% P₂O₅) was used as the phosphorus resource. For the latter, the full dosage was applied during sowing, and N was distributed by 20%, 40%, and 40% during sowing and during the first and second irrigations, respectively. Soil preparation and fertilization were the same for both growth cycles. No potassium fertilizer application was made because soils in this region have high potassium content, with average values of 3030 kg ha^-1 at a depth of 0.30 m [14].

The sowing was made by hand on dry soil on December 13, 2017 (Experiment 1) and on December 13, 2018 (Experiment 2) on a total surface of 2240 m². The seeding rate was 50 kg ha^-1 with germination of 85%. Each experimental plot was set in 12 rows, each 10 m long with a separation between lines of 18 cm (21.6 m²). The useful plot was 5 m long on the 10 central furrows (9 m²). Plants were subsequently thinned to leave a 160 m^-2 plant population density. On the same sowing date, a 150 mm irrigation depth was applied. During the cycle, three irrigation were made (37, 61, and 83 das), with an irrigation depth of 130 mm. The phenological stages when irrigation was applied after sowing were: rosette, elongation of the stem, and the initial formation of the capitulum. Weeds were controlled by hand with a hoe.

The yields for fresh forage and DM were determined at harvest. The content of DM was determined from a 0.4 m² sample randomly taken from the useful plot. For this purpose, a 0.74 m sample was taken from three of the central furrows of each plot. The sampled plants were weighted fresh, then pre-drying was performed for five days under the protection of a greenhouse, after which the samples were dried at 65 °C in a forced-air oven for 48 to 72 hours until a constant weight was attained [16, 17]. The DM yield was estimated by multiplying the fresh forage yield by the percent of DM forage in each useful plot. Leaf area index (LAI) was measured using AccuPAR model Lp-80 PAR/LAI Ceptometer (Decagon Devices, Inc., Pullman, WA, USA). The LAI measurements were collected before each one of the four harvest treatments. Three readings per plot were taken between 1200 and 1400 h solar time. Three measurements were made above and the other three below the canopy, parallel to the ground surface. The probe was positioned at a 45° angle with respect to rows. Plant height measured from ground level to plant apex of 10 random plants within each plot was determined before harvest. The plants sampled to estimate the DM content were used to analyze the nutritional value of the forage. The dry forage samples were ground in a Wiley® mill (Thomas Scientific, Swedesboro, NJ, USA) with a 1 mm mesh. The nitrogen content of each sample was determined using the Dumas Combustion Method number 990,03 of AOAC using Thermo Scientific Flash 2000 equipment, and the result was multiplied by 6.5 to obtain the percent of crude protein (CP) [18]. The Neutral Detergent Fiber (NDF) and the Acid Detergent Fiber (ADF) were analyzed in accordance with Goering and Van Soest [19]. The “in vitro” DM digestibility (IVDMD) was obtained from a ground dry forage sample placed in a Daisy incubator [20]. The content of net energy for lactation (NE_L) was obtained following the methodology proposed by Weiss et al. [21]. Yields of CP (CPY) and NE_L (NE_LY) per hectare were obtained by multiplying the contents of CP and NE_L by the DM yield per hectare estimated for each experimental plot. Digestible dry matter yield (DDMY) per hectare was obtained by multiplying IVDMD of forage by the DM yield per hectare.

Prevailing weather conditions for both years and the average values of 30 years during the same growing period are shown in Fig. (1). Meteorological data were obtained from a weather station located at the experimental site.

Fig. (1). Monthly temperature (a) and precipitation (b) during the safflower growing season in the two years of study and the average values of 30 years (1990-2019) at the La Laguna Experimental Station, Mexico.

The data obtained from this study were analyzed by using PROC MIXED from SAS version 9.3 [22]. A combined analysis of the data was performed using a randomized complete blocks design with four repetitions in a 4 × 4 factorial arrangement. Levene homogeneity of variances test was used to test the homogeneity of the variance, and it was found that the variances were homogeneous (P > 0.05). The analysis considered the block and the year as random variables, while the treatments were considered as fixed variables. Yields of DM, CP, NE_L, DDM, and concentrations of CP, NDF, ADF, and IVDMD were considered as fixed variables. The Tukey Kramer test was used at a level of P ≤ 0.05 to verify measurements. Also, a linear regression analysis was performed (P < 0.05) in order to determine the relationship between NDF and IVDMD with CP, NE_L, and IVDMD in the forage.

The combined analysis of the data using the year as a random effect shows no interaction between phenological stages and cultivars in all evaluated variables (P > 0.05) (Table 1). Then, the interaction was removed from the final model. The nutritional composition and the nutritional yield of safflower were affected by both the phenological stages and the cultivars significantly (P < 0.05) (Table 1).

The spineless cultivars Selkino and CD868 showed better or equal nutritional composition than the one observed for the spiny cultivars Gila and Guayalejo. All cultivars showed similar content of ADF and NDF. Regarding CP concentration, the better values occurred in Selkino and Guayalejo, while Selkino and Gila were outstanding in NE_L and IVDMD concentrations (Table 2). These differences between safflower cultivars were not observed in a study by Reta et al. [23], where a higher ADF concentration was reported in CD868 cultivar (372 g kg^-1) than that observed in the Selkino cultivar (341 g kg^-1).

Crude protein concentrations in all safflower cultivars (222 to 247 g kg^-1) were higher than those observed (119 g kg^-1 to 156 g kg^-1) in other studies, where spiny and spineless cultivars were evaluated [7, 12, 23-25]. Other chemical composition parameters of forage such as NDF (467 g kg^-1 to 477 g kg^-1), ADF (339 g kg^-1 to 390 g kg^-1), and NE_L (5.23 MJ kg^-1 to 5.73 MJ kg^-1 DM) were similar to those found in other studies, with values of NDF from 421 g kg^-1 to 462 g kg^-1, 331 g kg^-1 to 374 g kg^-1 of ADF [7, 12, 23-26] and from 5.15 MJ kg^-1 to 5.52 MJ kg^-1 DM of NE_L by Reta et al. [23]. Also, the levels of IVDMD observed in the current research (640 g kg^-1 to 674 g kg^-1) were similar to those observed by Landau et al. [7] and Ochoa-Espinoza et al. [24] in a spineless cultivar. Considering the nutrient composition of safflower, Landau et al. [7] (2004) indicated that forage of this crop could be used as silage for producing dairy cows.

Fig. (2). Relationship between Neutral Detergent Fiber (NDF) and In-Vitro Dry Matter Digestibility (IVDMD) with Protein Content (CP), Net Energy for Lactation (NEL), and IVDMD in the forage of four safflower cultivars harvested in four phenological stages (E) in the growing seasons 2017-2018 and 2018-2019.

The spineless cultivar CD868 produced DM yields and nutrients similar to those obtained from the spiny cultivars Gila and Guayalejo, while the spineless cultivar Selkino had smaller yields of DM, CP, and NE_L than the aforementioned three cultivars. However, Selkino produced a similar DDMY to the other three cultivars (Table 2). Differences in nutrient yields between Selkino and the other cultivars were associated with variations in agronomic traits. It is important to note that Selkino had a lower plant height at harvest (93 cm) compared with that achieved in the other cultivars (108 to 113 cm). Under this condition, Selkino was able to accumulate more DM in leaf (53%) and less in stem (47%) than the other three cultivars, where the accumulation of DM in the stem was 58 to 60% and in leaf from 40 to 42%. This behavior potentially decreased its DM yield but increased its IVDMD. Therefore, Selkino equaled the DDMY obtained by the other three cultivars. Landau et al. [7] reported a higher in vitro DM digestibility in leaves (729 g kg^-1) than in stems (546 g kg^-1) of safflower forage.

Fig. (3). Relationship between phenological stage (E), Leaf Area Index (LAI), and dry matter yield of four safflower cultivars harvested in four phenological stages in the growing seasons 2017-2018 and 2018-2019.

As the safflower harvest was carried out at a later phenological stage, the nutritional composition of the forage was less due to reductions in CP content and increases in fiber content. The parameter with a greater variation in the nutritional composition was CP (271 g kg^-1 to 194 g kg^-1), but its values were higher (152 g kg^-1 to 83 g kg^-1 DM) [23] or similar (272 g kg^-1vs. 125 g kg^-1 DM) than those reported by others studies [26, 27]. Regarding fiber concentration, the values for NDF increased (P < 0.05) from 418 g kg^-1 to 513 g kg^-1, while the values of ADF increased from 351g kg^-1 to 416 g kg^-1. This higher fiber content in the phenological stages with greater maturation was related to decreases in the content of CP, IVDMD (712 g kg^-1 to 591 g kg^-1), and NE_L (6.02 MJ kg^-1 to 4.94 MJ kg^-1 DM) in the forage (Fig. 2). In the last phenological stage evaluated in the study (E61), safflower forage was characterized by contents of CP of 194 g kg^-1, NDF of 513 g kg^-1, ADF of 416 g kg^-1, IVDMD of 591 g kg^-1, and 4.94 MJ kg^-1 of NE_L. The nutrient composition of safflower forage is considered acceptable when compared to that obtained in other traditional forages such as oat harvested at the beginning of heading in north-central México. At the heading stage, oat forage is characterized by containing 112 g kg^-1 CP, 349 g kg^-1 ADF, 547 g kg^-1 NDF and 5.4 MJ kg^-1 DM of NE_L [28].

Regarding this response of the safflower, Corleto et al. [27] indicated that the decrease in CP content and the increase in fiber concentration was associated with a lesser leaf proportion and a larger stem in the biomass of the crop in the harvests that were carried out in the phenological stages with greater maturation. Another study states that the increase in fiber in safflower forage in the phenological stages with greater maturation was due to the translocation of soluble cellular content from the leaves and stems to the seeds [26]. In other forages, the relation between the digestibility of the fiber and the nutritional composition of the forage showed a positive relationship between fiber digestibility and CP content [29, 30]. However, fiber digestibility has been negatively related to the concentrations of ADF and NDF in the forage [30, 31].

Yields of DM, CP, NE_L, and DDM presented a linear response with significant increases at every phenological stage when the harvest was delayed from E50 to E61. The greater yields for DM (10816 kg ha^-1), CP (2071 kg ha^-1), NE_L (52978 MJ ha^-1) and DDM (6350 kg ha^-1) occurred (P < 0.05) at the phenological stage E61 (Table 3). Delaying harvesting at stage E59 increased DM yield up to 8292 kg ha^-1, while the harvest at stage E61 allowed for DM yield up to 10816 kg ha^-1. The linear increase of DM yield when safflower grows from E50 (70 das) to E50 (76 das) is clearly associated with the higher DM accumulation ratio (265 kg ha^-1 day^-1) at the physiological stage of E61 by the crop. Flemmer et al. [15] indicated that the growth rate of safflower is related to stem elongation besides the formation and growth of branches and inflorescences after the E50 stage.

Furthermore, as the plant grows, the expansion of the leaf area also increases, which improves the interception of solar radiation. The leaf area increment in the present study is in line with the quadratic response of LAI observed between the stages E50 and E61 (Fig. 3), which also influenced the increase of DM yield. Reta et al. [23] observed a similar DM yield in safflower (9336 kg ha^-1) harvested at the E61 phenological stage. In contrast, Cazzato et al. [25] reported a higher DM yield (11600 kg ha^-1) than the DM production found at E61 when safflower was harvested at the 25% flowering stage. The increase of nutrient yields when harvesting from E50 to E61 follows the same pattern as the increase of DM yield when harvesting at the same physiological stages. Therefore, the yields of CP, NE_L, and digestible DM are clearly explained by the effect of DM yield since the percentage of these nutrients in forage decreased as the harvest progressed from E50 to E61. This response allowed safflower to produce the highest nutrient yields at the E61 stage.

The null difference in the interaction cultivars × harvest age indicates that the greater yields of DM obtained in the study may be attained with spiny and spineless cultivars. However, when spiny cultivars are used, the highest forage yield free of spines is obtained during the phenological stage E50, since in later stages, the spines on the leaves and inflorescences are already present. The forage in phenological stage E55 may be considered acceptable since the spines present on the plant are not yet fully developed. This is important because due to the spiny nature of traditional safflower cultivars, farmers can decide whether or not to use this forage in the diet of dairy cattle since cows are more susceptible than goats and sheep to ulcerations of the mouth caused by spines present in traditional safflower forage [32].

The use of spineless safflowers allowed the production of forage free of spines until the last phenological stage (E61) evaluated in the study. Since the nutritional composition of the safflower forage maintained an acceptable level during the two last phenological stages (E59 and E61), the increase in DM yield at phenological stage E61 (58%) allowed the spineless safflowers to increase their yields by 29% for CP, by 39% for NE_L, and by 41% for digestible DM as compared to the yields of spiny safflower harvested during phenological stage E55, which is the latest stage when forage with undeveloped spines may be obtained.