text:
Novel Insight into the Volatile Profile and Antioxidant Properties of Crocus sativus L. Flowers

The current production system of saffron spice generates hundreds of tons of waste. Thus, the aim of this study was to value both saffron and its floral by-products as a source of natural bioactive extracts, studying the in vitro antioxidant capacity, the composition of the volatile fraction by GC-MS/MS, and the determination of crocetins esters by HPLC-PDA. Saffron stigmas and floral by-products showed a high content of polyphenols and different antioxidant properties. Floral bio-residues (tepals, stamens, and styles) presented a high concentration of anthocyanins, and stigmas had high levels of flavonoids, β-carotene, and total crocins. In stigmas, 25 different volatile components were found, with safranal the most relevant. Floral by-products volatile composition consisted of 55 compounds with varying amounts depending on the drying treatment; all the samples presented acetic acid, 2(5H)-furanone, and phenylethyl alcohol. Therefore, saffron stigmas and flower by-products represent a sustainable source of bioactive ingredients for innovative healthy food formulations.

Introduction
Saffron (Crocus sativus L.), a traditional Mediterranean plant, is a perennial herb that belongs to the Iridaceae family. It is employed as a spice, natural colorant in food, and a flavoring agent. This spice is obtained from the flower of C. sativus, which is composed of three golden yellow stamens, six purple tepals, and one red pistil that culminates with three red branched stigmas, whose length exceeds that of the tepals, which, when dried up, constitute the saffron spice [1]. Only flower stigmas are used for saffron production, while tepals and stamens are simply discarded. For the production of 1 kg of saffron, around 350 kg of tepals are generated as a by-product. Thus, the current production system is not sustainable since it generates hundreds of tons of waste, with a high environmental impact [2].
In addition to its organoleptic properties, saffron is occasionally used as a textile dye, in perfumes, and in medicine due to its therapeutic properties since it has been considered for centuries as a medicinal plant in many cultures [3]. The therapeutic activity of saffron is mainly due to its valuable bioactive compounds (carotenoids, terpenes, and flavonoids). The main chemical components are crocin, crocetin and safranal. Crocetin, the principal apocarotenoid in saffron, is the precursor of crocin, and the most abundant compound is trans-crocetin di(β-D-gentiobiosyl) ester. The volatile compounds, including terpenoids, along with phenolic compounds that are widely distributed in plants, are responsible for their sensory properties. In saffron, 3-Cyclohexadiene-1-carboxaldehyde, 2,6,6-trimethyl-, known as safranal, is the major volatile compound contributing to saffron aroma [4].
One of the greatest interests of these bioactive compounds is due to their high antioxidant capacity and free radical scavenging activity. The generation of reactive oxygen species is a normal process in cells, but the uncontrolled generation and concomitant increase in reactive oxygen species (ROS) level in the body results in "oxidative stress", which is considered to be the main cause of various diseases [5]. The antioxidant substances play a role in protecting biological systems against the effects of oxidative processes on macromolecules. Many of those substances are plant-derived natural molecules, such as carotenoids or polyphenolic compounds, that could contribute to preventing and fighting against several diseases in which ROS are involved [6].
Currently, some studies have reported that tepals, considered waste, contain valuable bioactive compounds such as flavonoids and polyphenols with potential antioxidant activity [7]. Therefore, saffron floral by-products could be natural antioxidant sources to be used in food, increasing the resource efficiency and the competitiveness of this sector. This waste raw material could be valorized through the development of innovative high added-value food ingredients, and this could potentially increase saffron demand on the market.
The main objective of this study was to value both saffron and its floral by-products through their use as a source of bioactive extracts, studying the in vitro antioxidant capacity, the metal chelating activity, the composition of volatile compounds by GC-MS/MS, and the determination of crocetins esters and total crocins by HPLC in saffron and its floral byproducts. As far as we know, no research has been published to date that delves into the composition of the volatile fraction in saffron floral by-products and into the changes in the volatile profile depending on the drying treatment used. Moreover, the results of this research would support the production of saffron as a more sustainable and profitable agronomic resource.
For all the experiments, ultrapure Milli-Q water (Millipore Corp., Bedford, MA, USA) was used.

Plant Material
Saffron floral by-products (SFL1, SFL2) were obtained from different producers in Castilla-La Mancha region (Spain) during the 2020 harvest season and cultivated following the requirements established by the Protected Designation of Origin "La Mancha Saffron" according to DOCM [8]. The flowers were composed of all the parts of saffron flowers (tepals, stamens, and styles), except for the stigmas that were detached manually from the rest of the flower after being harvested, following the traditional procedure.
All fresh flowers were frozen in liquid nitrogen and stored at −80 • C until freeze-dried in a freeze-dryer Christ Alpha 2-4 (B. Braun Biotech International, Melsungen; Germany) for 48 h to constant weight (initial temperature −25 ± 2 • C and pressure 0.220 mbar). Then, they were crushed and sieved through a 500 µm mesh size and stored at −20 • C until further analysis. For the volatile composition analysis, in addition to fresh and freezedried samples, oven-dried samples were studied, which were dried for 24 h at 60 • C in an air oven.
Saffron stigmas were supplied dehydrated by the Spanish company Verdú Cantó Saffron Spain, S.L and were from Spanish (Castilla-La Mancha region), Greek (Kozani area), and Iranian cultivations (Torbat zone). The moisture of all samples was lower than 11%, and saffron threads were crushed and sieved through a 500 µm mesh size and stored at 4 • C until further analysis.

Polyphenolic Compounds
The extractions were prepared using a methanol solution and a sample/methanol ratio 1:20 (w/v) for the freeze-dried flowers and a ratio 1:50 (w/v) for saffron stigmas. The extracts were shaken for 1 h in the dark at 400 rpm with a magnetic stirrer (Ovan, mod. MultiMix Heat D-MMH30E) and then sonicated for 15 min and centrifuged at 10,000× g for 10 min at 4 • C. Then, the supernatants were filtered (0.45 µm PTFE filter, Millipore, Madrid, Spain) and stored at −20 • C. All extractions were performed in triplicate.

Beta-Carotene
The extractions were prepared using n-hexane:acetone:ethanol (50:25:25, v/v/v) solution and a sample/solvent ratio 1:20 (w/v) for the freeze-dried flowers and a ratio 1:50 (w/v) for saffron stigmas. The extracts were shaken for 10 min at 400 rpm with a magnetic stirrer (Ovan, mod. MultiMix Heat D-MMH30E), keeping them in ice, and then centrifugated at 14,460× g for 20 min at 4 • C. Then, the supernatants were filtered (0.45 µm PTFE filter, Millipore, Madrid, Spain) and stored at −20 • C. All extractions were performed in triplicate.

Anthocyanins
The extractions were prepared using 50% ethanol (0.1% HCl) solution and a sample/solvent ratio 1:20 (w/v) for the freeze-dried flowers. The extracts were sonicated for 30 min and centrifugated at 10,000× g for 15 min. Supernatants were filtered (0.45 µm PTFE filter, Millipore, Madrid, Spain). All extractions were performed in triplicate.

In Vitro Antioxidant Properties
For the study of the antioxidant properties, the extracts obtained from polyphenolic compounds extraction (Section 2.3.1.) were used.

2,2-Diphenyl-1-picrylhydrazyl (DPPH) Free Radical Scavenging Method
The free radical scavenging activity was determined using DPPH method, following the methodology from Brand-Williams et al. [9], with some modifications. The DPPH radical was prepared dissolving 0.0035 g with 10 mL of methanol. The mixture was shaken and kept in dark for 30 min. The absorbance decrease was measured at 515 nm (UV/Vis Spectrophotometer T80; PG Instruments Limited, UK). Trolox (10 mM) was used as a reference standard at different concentrations (0.50-4.00 mmol/L). The results were expressed as mmol of Trolox Equivalents (TE) per 100 g dw (dry weight) of sample.

2,2 -Azinobis (3-ethylbenzothiazoline-6-sulfonic acid) (ABTS) Free Radical Scavenging Method
The ABTS cation radical method for measuring the antioxidant capacity was performed by adapting the methodology from Re et al. [10]. ABTS radical was prepared mixing ABTS (7 mM) with K 2 S 2 O 8 (2.45 mM), and they reacted for 16 h in dark at room temperature. Subsequently, the solution was diluted with ultrapure water until its absorbance was adjusted to 0.70 ± 0.02 at 734 nm. Trolox (10 mM) was used as a reference standard in different concentrations (0.20-3.00 mmol/L). The results were expressed as mmol of TE per 100 g dw of sample.

Ferric Reducing Antioxidant Power (FRAP) Method
The FRAP method adopted from Benzie and Strain [11] was used. Briefly, the FRAP reagent was prepared fresh daily by mixing 300 mmol/L acetate buffer (pH 3.6), 10 mmol/L TPTZ solution in 40 mmol/L HCl, and 20 mmol/L FeCl 3 ·6H 2 O solution in a volume ratio of 10:1:1, respectively. The absorbance was measured at 593 nm. Trolox (10 mM) was used as standard solution in the range of 0.01−5.00 mmol/L. The results were expressed as mmol of TE per 100 g dw of sample.

Total Polyphenols Content (TPC)
The total polyphenols were determined using the Folin-Ciocalteu methodology following the method described by Singleton et al. [12]. The Folin-Ciocalteu reagent was mixed with ultrapure water 1:10 (v/v). Gallic acid (1 mM) was used as a reference standard in the range of 0.00-4.72 mg/L. For the assays, 100 µL of the different extracts were mixed with 400 µL of phosphate buffer (50 mM) at pH 7.8, and 2.5 mL of Folin-Ciocalteu reagent was added. After 2 min, 2 mL of Na 2 CO 3 (75 g/L) were added and kept at 50 • C for 10 min. The absorbance was measured at 760 nm in a spectrophotometer (UV/Vis), and the results were expressed as mg Gallic Acid Equivalents (GAE) per g dw of sample.

Total Flavonoids Content (TFC)
The total flavonoid content was determined as described by Çam and Hışıl [13]. Catechin was used as standard; in order to calculate the flavonoid content, a calibration curve ranging from 20-100 mg/L was prepared. The absorbance was measured at 510 nm, and the results were expressed in mg of Catechin Equivalents (CE) per g dw of sample.

Total Anthocyanins Content (TAC)
The total anthocyanin content was determined by the method from Figueira et al. [14] with some modifications. Delphinidin was used as an external standard, and the total anthocyanin content was calculated using a calibration curve (10-150 mg/L). The final absorbance was determined by the difference between the measured reference absorbance and the measured sample absorbance. The absorbance was measured at 520 nm, and the results were expressed in mg of Delphinidin Equivalents (mg DE) per g dw of sample.

Total Beta-Carotene Content
The change in absorbance was measured using a UV-Vis spectrophotometer at 450 nm. Carotene content was determined using the Beer-Lambert law. The carotene concentration was calculated using the extinction coefficient of beta-carotene in hexane (2.505 M −1 cm −1 ). The results were expressed as mg of beta-carotene per 100 g dw of sample (%).

Iron (II) Chelating Activity
For the metal chelating activity, the extractions were performed in the same procedure previously explained for the extraction of polyphenolic compounds but using ultrapure water as extracting agent.
The chelation of iron (II) ions was performed as described by Carter [15]. Briefly, a final concentration of 0.02 g/mL of extract were used and mixed with 100 µL of 2.0 mM aqueous FeCl 2 and 900 µL methanol. After an incubation of 5 min, the reaction was initiated adding 400 µL of 5.0 mM ferrozine. Then, after a 10 min equilibrium period, the absorbance was measured at 562 nm. The iron chelation activities were calculated from the absorbance of the control (Ac) and of the sample (As) using the following equation and expressed as % inhibition:

Color
The color was measured with a Minolta CR-300 Chroma Meter (Japan) colorimeter, using L*, a*, b* scale (CIELAB system). The results were expressed as luminosity L*, a* (greenness/redness), b* (blueness/yellowness), Hue angle (h) and Chroma (C*), which were calculated according to the following equations, respectively:

Extract Preparation
The aqueous extracts of saffron stigmas and flowers were prepared according to ISO 3632 [16]. Briefly, 0.5 g of sample was mixed with 900 mL of ultrapure water in a 1 L volumetric flask. The solutions were shaken for 1 h in the dark at 1000 rpm on a magnetic stirrer. The flask was filled to 1 L and homogenized through agitation. Then, 20 mL of solution were transferred to a 200 mL volumetric flask which was filled with ultrapure water. The extractions were filtered (0.45 µm PTFE filter, Millipore, Madrid, Spain) and transferred into a vial for HPLC analysis. All extractions were performed in triplicate.

HPLC-PDA
The identification and quantification of crocetin esters were carried out by highperformance liquid chromatography (HPLC) using HPLC Altus TM 10 PerkinElmer (Waltham, MA, USA) equipped with a C18 KromaPhase column (150 × 4.6 mm inner diameter, 3.5 µm) (Scharlab, Barcelona, Spain) that was equilibrated at 30 • C. The analysis were carried out as described by Valle García-Rodríguez et al. [17] with some modifications. The eluents used were acetonitrile (A) and water (B), with a proportion of 20% A and 80% B. The flow rate was 0.8 mL/min, and the injection volume of saffron extracts was 25 µL. Crocetin esters were measured at 440 nm in UV-VIS (UV-6300PC double beam spectrophotometer) with the Photodiode Array (PDA) Detector (PerkinElmer, Waltham, MA, USA).

Volatile Composition
The volatile composition was determined using headspace solid phase micro-extraction (HS-SPME). After several preliminary tests to optimize the extraction system, accurately amounts between 20 and 40 mg of saffron stigmas and between 120 and 140 mg of saffron flowers were weighted and added into a 40 mL vial with polypropylene caps and PTFE/silicone septa, and isoamyl acetate (1000 mg/L, internal standard for semi-quantification of compounds). Fresh, freeze-dried, and air oven-dried samples of saffron flowers were used in the analyses. Then, the vial was placed in an AOC-6000 Plus autosampler (Shimadzu Corporation, Kyoto, Japan), and after 5 min of equilibration time, a 50/30 µm DVB/CAR/PDMS fiber (1 cm) was exposed to the sample headspace for 45 min at 40 • C (with agitation, 250 rpm).
The separation and identification of compounds were performed by GC2030 (Shimadzu Scientific Instruments, Inc., Columbia, MD, USA) in a Sapiens X5MS column (Teknokroma, Barcelona, Spain), 30 m × 0.25 mm i.d., 0.25 µm f.t., and coupled with a mass spectrometer detector (TQ8040 NX triple quadrupole mass spectrometer; Shimadzu Scientific Instruments, Inc., Columbia, MD, USA). Only the single quadrupole acquisition mode was used on the TQ8040 NX (Q3 Scan; scan speed 5000 amu/s; mass range 40-400 m/z; event time 0.100 s). The oven temperature program was as follows: (i) initial temperature of 35 • C and hold 5 min; (ii), increment of 5 • C/min up to 150 • C/min, and hold 1 min; (iii) increment of 10 • C/min up to 280 • C and hold for 5 min. Helium column head pressure was 47.6 kPa (constant linear velocity mode of 36 cm/s). Injector, ion source, and interface were at 250, 230, and 280 • C, respectively. Helium was used as gas carrier, column flow 1 mL/min, with split ratio 1:50 and purge flow of 6 mL/min.
Retention indexes of a commercial alkane standard mixture were used to identify the compounds, as well as the National Institute of Standards and Technology (NIST) 17 Mass Spectral and Retention Index Libraries. The identification was considered tentative when it was based only on mass spectral data, and only compounds with spectra similarity > 90% were considered correct hits. Linear retention similarity filter was set at ±10 units. This volatile compound extraction method has been previously used for the analysis of different food matrices, according to Clemente-Villalba et al. [18].

Statistical Analysis
All determinations were performed in triplicate. Results were expressed as the mean ± standard deviation. The mean comparisons were carried out using an analysis of variance (ANOVA) and by the Tukey multiple range test using SPSS version 21.0 software package (SPSS Inc., Chicago, IL, USA). The significant differences were established as p < 0.05.

Antioxidant Properties, Bioactive Content and Iron (II) Chelating Activity
Polyphenols, which are plant-derived natural molecules found as secondary metabolites, present important biological activities, with their antioxidant capacity one of the most important properties for physiological function. The total polyphenols and flavonoid content, as well as the antioxidant effects, were studied in saffron and its floral by-products extracts.
The results of the antioxidant activity and bioactive content are shown in Table 1. Saffron floral by-products and saffron stigmas presented a high concentration of total polyphenols, between 32-36 mg GAE/g dw, except for Spanish saffron, which had a significant highest amount (44.80 ± 2.30 mg GAE/g dw). Regarding the total flavonoid content, saffron stigmas from Spain, Iran, and Greece showed high levels of flavonoids in the range of 15-18 mg CE/g dw. However, saffron flowers had a lower concentration of total flavonoids (4-5 mg CE/g dw) than saffron stigmas.
These results obtained for saffron flowers were in accordance with the ones described by Sun et al. [19], reporting values about 30 mg/g for TPC and values lower than 10 mg/g for TFC in saffron tepal methanol extracts from China. Moreover, the values obtained for saffron stigmas were higher than those reported in the study of Karimi et al. [6], in which saffron stigmas methanol extracts from Iran presented 6.5 ± 0.02 mg GAE/g dw for TPC and 5.8 ± 0.12 mg rutin equivalents/g dw for TFC. 1 Means ± standard deviation in the same line followed by different lowercase letters indicate statistically significant differences at p ≤ 0.05 for each sample (n = 3). 2 The antioxidant capacity is expressed as mmol Trolox Equivalents (TE) per 100 g dw of sample. 3 The Total Polyphenols Content (TPC) are expressed as mg Gallic Acid Equivalents (GAE) per g dw of sample. 4 The Total Flavonoids Content (TFC) are expressed as mg of Catechin Equivalents (CE) per g dw of sample. 5 The Total Anthocyanin Content (TAC) were expressed in mg of Delphinidin Equivalents (mg DE) per g dw of sample. 6 The total bta-carotene is expressed as mg of beta-carotene per 100 g dw of sample (%). SFL1, SFL2: Saffron floral by-products from two different producers; n.d.: not determined.
Within the group of flavonoids, anthocyanins are one class of pigments that also present antioxidant properties [20]. Saffron flowers without stigmas, SFL1 and SFL2, contained high levels of anthocyanins, showing SFL2 a significantly higher concentration than SFL1. Serrano-Díaz et al. [21] also found the highest anthocyanin content in saffron tepals from Spain and were also detected in whole flowers and floral bio-residues.
Regarding carotenoids, they are organic pigments from the group of isoprenoids that are found naturally in plants and exert several beneficial functions in the body due to their antioxidant properties, among others. Based on their structures, carotenes are one of the main subclasses of carotenoids, with β-carotene and α-carotene the two major types. Beta-carotene was found to be the precursor of crocins. Crocins, which are glycosyl esters of crocetin, are the main chemical components of saffron, responsible for many of its pharmacological and biomedical properties. The synthesis of these apocarotenoids involves several reactions, including the cleavage of beta-carotene and zeaxanthin, oxidation and glycosylation steps [22].
Saffron stigmas showed a high proportion of total beta-carotene being higher than 70%. Greek saffron presented the highest concentration (90.82 ± 8.37%) and Iranian saffron the lowest, having 20% less total beta-carotene than Greek saffron (71.06 ± 6.47%) ( Table 1). With respect to the crocetin esters content, the trans isomers, which are the majority, were studied. Trans-4-GG crocin was found in higher concentrations than trans-3-Gg-crocin in all saffron samples. Spanish stigmas presented the significant highest proportion of total crocins (17.04 ± 0.18%) and of trans-4-GG crocin (11.91 ± 0.13%) ( Table 2). The same tendency for beta-carotene content was shown by Iranian stigmas, having the lowest proportion of total crocins (13.19 ± 0.78%) and of trans-4-GG crocin (9.22 ± 0.30%) ( Table 2). These results were in accordance with those obtained by Valle García-Rodríguez et al. [17], that have reported the quantification of these two crocins (trans-4-GG and trans-3-Gg) in saffron from Italy, Iran, Greece, and Spain by HPLC. Furthermore, the values of Spanish saffron were consistent with the content of crocetin esters reported by Moratalla-López et al. [23], in which these compounds represented 16-28% of the saffron composition. Thus, the level of apocarotenoids may vary because of different geographical origins, processing, and storage conditions. Saffron flowers presented low amounts of total beta-carotene (28-39%) that are mainly located in the yellow stamens. Nevertheless, beta-carotene was reported in higher quantity in other Crocus species such as C. ancyrensis also due to the yellow color of the tepals [22] ( Table 1). Moreover, in saffron floral by-products, crocins were only detected in SFL1 but in very low amounts (0.042%) due to the presence of small fragments of stigmas remaining after their detachment (Table 2). However, Rubio Moraga et al. [22] reported that crocins in C. sativus saffron floral by-products from Spain were not detected.
DPPH, ABTS, and FRAP assays results revealed the different antioxidant activities of saffron and its floral by-products (Table 1). DPPH and ABTS tests evaluate the in vitro antiradical activity, and FRAP assay the reducing potential of the extracts. Saffron stigmas presented stronger antioxidant capacity than the floral by-products, which could be related to the higher amounts of total flavonoids and crocins content since previous studies have demonstrated the significant antioxidant activity of flavonoids and crocins [24]. Nevertheless, saffron flowers, SFL1 and SFL2, showed a good antioxidant ability (by DPPH, ABTS, and FRAP) that could be related to their high polyphenols content. These results are similar to other studies that reported the antioxidant capacity of Indian and Iranian saffron stigmas [6,25], the potential antioxidant activity of commercial saffron powder from the large-scale Italian market [26], and also to that of Sun et al. [19] in which it is revealed the strong antioxidant ability of saffron tepals from China.
Regarding the ability of the saffron extracts to chelate iron (II), the results are shown in Figure 1. Bivalent transition metal ions, such as iron, play an important role as catalysts of oxidative processes, participating in hydroxyl radical generation via the Fenton reaction [27]. Thus, excess metal ions could lead to the formation of free radicals generating high levels of oxidative stress, but these processes can be delayed by iron chelation. Apart from the ferric reducing power activity (FRAP) assay determination that indicated the reduction potential of Fe 3+ to Fe 2+ , the iron (II) chelating activity of saffron stigmas and flowers was also studied.  Therefore, the chelating activity of Fe 2+ and antioxidant capacity that presented saffron and its floral by-products could prevent oxidative damage, protecting against oxidative stress. These abilities may be linked to their bioactive compounds content that were Saffron stigmas presented a good chelating activity, around 60% of inhibition, with the highest for the Greek saffron ( Figure 1). This ability could be due to their flavonoid content since previous research demonstrated that these bioactive compounds have the capacity to neutralize reactive radicals and sequestrate metal ions, suppressing Fenton reactions [28]. Saffron flowers showed a low iron chelating activity (less than 20% of inhibition), with the capacity of SFL2 significantly higher than SFL1 (Figure 1). This fact might be related to the anthocyanins content of SFL2 (Table 1). These results were similar to those reported by Sánchez-Vioque et al. [29], indicating that saffron tepals from Spain had very low chelating activity.
Therefore, the chelating activity of Fe 2+ and antioxidant capacity that presented saffron and its floral by-products could prevent oxidative damage, protecting against oxidative stress. These abilities may be linked to their bioactive compounds content that were located differently on the flower; stigmas had carotenoids and flavonoids, and the rest of the flower presented anthocyanins and other polyphenolic compounds.

Color
The color parameters of saffron and its floral by-products are shown in Table 3. The L* values were around 30 in all samples indicating a low luminosity and lightness. Regarding powdered saffron floral by-products, SFL1 and SFL2, positive a* values and negative b* values represented a color in the ranges of red and blue, respectively, which were characteristic of the violet color of tepals. It should be noted that the b* value of SFL2 (−4.77 ± 0.14) indicated a higher blue color intensity with respect to SFL1 (−1.63 ± 0.25). This fact may be mainly related to the high concentration of anthocyanins pigments found in SFL2, which render a blue, red, or purple color [30]. Moreover, h values ( • ) of SFL1 (334.50 ± 2.46) and SFL2 (313.00 ± 0.19) represented a tone in the blue/purple range, and the color was more saturated and intense in SFL2 than SFL1, represented by the chroma (C*) values (6.53 ± 0.18 and 3.78 ± 0.22, respectively). Saffron stigmas from Spain, Iran, and Greece showed positive a* and b* values, representing colors ranging from red to yellow, respectively, which are characteristic of saffron spice. The h values were around 26-28 • which indicated a red-orange tone that was saturated and intense in all saffron samples (C* values around [10][11][12]. This yellow to red color observed in the stigmas of Crocus sativus L. was due to the presence of a high amount of carotenoids, such as crocetin and crocins, which are responsible for the coloring power of the saffron spice exhibiting red, orange, and yellow colors [22].
The results showed that the color parameters variability within the flowers and stigmas of Crocus sativus L. is related to the number of bioactive compounds with coloring properties in the different parts of the plant.

Volatile Composition
The evaluation of volatiles is an important aspect that contributes to the aroma of foods. Furthermore, volatile secondary metabolites also present important bioactive properties, such as antioxidant, antimicrobial, and anticancer activities [31].
The results of the volatile composition of saffron stigmas from Spain, Iran, and Greece are shown in Table 4. In all saffron samples, a total of 25 different volatile compounds were identified and quantified, with the more representative chemical families, ketones (n = 8), terpenes (n = 8), aldehydes (n = 3), esters (n = 3), alcohols (n = 2) and acids (n = 1). Regarding the concentration of the different volatile compounds, there were significant differences between the samples, except for safranal (1,3-Cyclohexadiene-1-carboxaldehyde, 2,6,6trimethyl-, V16) which was the compound found in a higher amount in all the samples regardless of the saffron origin (7443, 7656 and 7429 µg/g for Spanish, Iranian and Greek saffron, respectively). The concentrations of safranal were higher than those reported by other studies. Kanakis et al. [32] indicated maximum values of safranal in Greek saffron around 6879 µg/g, Culleré et al. [33] found 1365 ± 81.9 µg/g in Spanish saffron, and Jalali-Heravi et al. [34] showed values of 4356 µg/g for Iranian samples. These differences could be due to the post-harvest processing, such as the dehydration procedures since safranal is formed at elevated temperatures from picrocrocin (50-55 • C), but after this process, safranal could still be generated via trans-crocetin esters from other volatile compounds [4].
In summary, saffron from Greece had the highest concentration of total volatile compounds (around 13,067 µg/g) with respect to the Iranian and Spanish saffron samples (around 12,503 µg/g and 11,104 µg/g, respectively). These differences could be related to the geographical origin, which is considered an essential factor in the concentration of volatiles, among others, such as the harvest season, dehydration temperature and duration, and the storage time and conditions [39].
The aroma of saffron is mainly due to the contribution of aldehydes of saffron, such as safranal and its derivative 4-hydroxy-2,6,6-trimethyl-1-cyclohexene 1-carboxaldehyde, but other compounds could contribute as well, with the characteristic saffron aroma developed during the post-harvest treatment. With respect to the odor descriptors, the saffron aroma is especially herbal, spicy, sweet, fresh, woody, floral, and musty, among others (Table 4).
Regarding the flower by-products samples from Spain fresh and dried by two different methods (freeze-dried and air oven-dried), 55 volatile compounds were isolated, identified, and quantified ( Table 5). The identification parameter of these compounds, namely Kovats index, retention time, chemical family, and odor descriptors, were indicated in Supplementary Table S1.
The chemical families more representative found in all samples were: esters (n = 13), aldehydes (n = 12), and acids and alcohols (n = 5). All floral by-product samples presented acetic acid (V1), 2(5H)-furanone (V14), and phenylethyl alcohol (V31). SFL2 (air-oven dried) had the highest concentration of acetic acid (161.55 µg/g), followed by fresh SFL1 (145.21 µg/g). The production of this organic acid could occur by the conversion of some carbohydrates in the presence of oxygen.
SFL2 (air-oven dried) presented a significant highest level of 2(5H)-furanone (V14) (452 µg/g), while SFL1 (air-oven dried) showed the lowest concentration (4.11 µg/g). However, the generation of furanone, such as 2(5H)-furanone, may be due to several processes: via oxidation reactions, by microorganisms, or spontaneous formation via the Maillard reaction between sugars and amino acids during heating [40]. Thus, the content of this volatile compound could be related mainly to the sugar and amino acid composition of the fresh flowers SFL1 and SFL2, leading to different concentrations of 2(5H)-furanone in the air-oven dried flowers. Moreover, 2(5H)-furanone is authorized as a flavoring substance by the Joint FAO/WHO Expert Committee on Food Additives (JECFA) and is defined in its specification as a "rich winey meat-like aroma". This fact could be interesting for the use of dried saffron flowers as food ingredients. This compound has also been studied as bactericidal [41].
Phenylethyl alcohol was present in the SFL2 (air-oven dried) sample in high concentrations (35.63 µg/g). However, minimal concentrations were found in the SFL1 (air-oven dried) sample (0.49 µg/g). This volatile was also found in other flower extracts (rose, hyacinth, geranium), presenting a rose-like odor and antimicrobial and antifungal properties, and it is widely used in foods and cosmetics [42].
It should be noted that safranal (V39) was only found in SFL2 (fresh and air-oven dried) but in very low amounts (0.11 and 7.53 µg/g, respectively). This could be due to the presence of traces of saffron stigmas in the floral by-products, with the amount of safranal higher in SFL2 (air-oven dried) since it is generated at the temperatures used (50-55 • C) for drying. Moreover, isophorone (V35) was also present in SFL2 (22.28 µg/g) which may be produced from safranal and/or from other precursors formed in the oven during the heating process.
Other compounds that should be highlighted due to their high concentration were 3-hydroxy 2-butanone (acetoin, V4) and butanoic acid, 3-methyl (isovaleric acid, V9). These volatiles were found in fresh SFL1 at 299.42 µg/g and 177.52 µg/g, respectively. Acetoins, which exist widely in nature, are mainly found in higher plants that have the ability to synthesize them using different enzymes, but their biological mechanisms remain unclear; however previous studies have shown their antimicrobial actions [43]. Isovaleric acid is an isomer of valeric acid, which is a compound naturally present in plants as a metabolite that could be formed by the secondary metabolism of plants [44]. Moreover, these volatiles were only found in fresh samples, meaning that they may be degraded by heat.
In summary, some volatiles present in fresh flowers increased by air drying in the oven due to the thermal treatment, and new volatiles also appears in the air-dried flowers. Some volatile compounds (V3, V17, V26, V32, V44) were only found in SFL2 (air-oven dried) samples, which showed the highest concentration of total volatile compounds (833.49 µg/g). These results are justified by the generation of volatiles from other non-volatile precursors by heat processing. Thus, the post-harvest procedures, such as drying, are very important when considering the volatile quality in terms of the contribution of the volatile compounds that are generated during this processing step.

Conclusions
The studied Spanish, Iranian, and Greek saffron stigmas and Spanish saffron floral by-products could have a great potential to develop new high-added-value ingredients due to their antioxidant properties and bioactive content. Floral by-products showed a good concentration of total anthocyanins, and saffron stigmas had a high concentration of total flavonoids, β-carotene, and crocetin esters. Furthermore, among the volatile composition, the saffron studied presented high levels of safranal, which contributes to the bioactivity and aroma of this spice. The floral by-products volatile composition from fresh to dried samples was highly influenced by the drying method, of which some may present biological activities, especially antioxidant and antimicrobial properties. Therefore, these results suggest that saffron and its floral by-products are natural sources of antioxidant compounds that could be used as sustainable, innovative ingredients to apply in food for the development of novel functional food products or for other human health applications.
Supplementary Materials: The following supporting information can be downloaded at: https: //www.mdpi.com/article/10.3390/antiox11091650/s1, Table S1: Identification, quantification and descriptors of volatile compounds found in saffron floral by-products. Funding: This article is supported by the Partnership for Research and Innovation in the Mediterranean Area (PRIMA) program, project 174-SAFFROMFOOD. The PRIMA program is supported by the European Union. The research was also supported by the regional agency "Diputación de Alicante", Alicante, Spain. Grant EQC2018-004170-P funded by MCIN/AEI/10.13039/501100011033 and by ERDF A way of making Europe for GC-MS analysis.


provenance: AgriculturalAndFoodSciences,2022-2022-Agricultural And Food Sciences-2022.jsonl.gz.gz

text:
Construction of Chromosome Segment Substitution Lines and Inheritance of Seed-Pod Characteristics in Wild Soybean

Genetic populations provide the basis for genetic and genomic research, and chromosome segment substitution lines (CSSLs) are a powerful tool for the fine mapping of quantitative traits, new gene mining, and marker-assisted breeding. In this study, 213 CSSLs were obtained by self-crossing, backcrossing, and marker-assisted selection between cultivated soybean (Glycine max [L.] Merr.) variety Suinong14 (SN14) and wild soybean (Glycine soja Sieb. et Zucc.) ZYD00006. The genomes of these 213 CSSLs were resequenced and 580,524 single-nucleotide polymorphism markers were obtained, which were divided into 3,780 bin markers. The seed-pod-related traits were analyzed by quantitative trait locus (QTL) mapping using CSSLs. A total of 170 QTLs were detected, and 32 QTLs were detected stably for more than 2 years. Through epistasis analysis, 955 pairs of epistasis QTLs related to seed-pod traits were obtained. Furthermore, the hundred-seed weight QTL was finely mapped to the region of 64.4 Kb on chromosome 12, and Glyma.12G088900 was identified as a candidate gene. Taken together, a set of wild soybean CSSLs was constructed and upgraded by a resequencing technique. The seed-pod-related traits were studied by bin markers, and a candidate gene for the hundred-seed weight was finely mapped. Our results have revealed the CSSLs can be an effective tool for QTL mapping, epistatic effect analysis, and gene cloning.

INTRODUCTION
Cultivated soybean (Glycine max [L.] Merr.) is an important food crop and oil crop in the world, and it is the most important source of plant protein (Nadeem et al., 2021). Therefore, the development of soybean industry is of great significance. Good variety is an important prerequisite for high yield and high yield of soybean. According to various traits of soybean, cultivating good varieties and improving yield and quality of soybean are the goals of breeding workers (Reinprecht et al., 2006). The emergence of molecular marker-assisted selection (MAS) and QTL mapping has significantly shortened the breeding time (Hasan et al., 2021). Some June 2022 | Volume 13 | Article 869455 Zheng et al.
CSSLs Construction and Inheritance Study excellent loci and genes that are of great significance to breeding have been identified by populations, such as F 2 , BC 1 , DH, RIL, and IF 2 (Mian et al., 1996;Song et al., 2005;Reinprecht et al., 2006;Alfred et al., 2012;Botero-Ramirez et al., 2020). The genetic composition of each member in the above population is diverse, which can provide abundant genetic variation information and calculate multiple genetic effects. However, the diversity of genetic background will also bring some interference to QTL mapping and affect its accuracy Qiao et al., 2016). The CSSLs constructed based on backcrosses have similar genetic backgrounds, which can eliminate interference from genetic noise and improve the accuracy of QTL mapping (Ali et al., 2010;Balakrishnan et al., 2019). CSSLs are constructed by continuous backcross between recurrent parent and donor parent. The genetic background is clean and uniform, and only a few or a single donor genome introgressed fragments are different among the lines. CSSL population can accurately locate QTL loci for target traits. In addition, the process of backcrossing can introduce valuable genes into cultivated varieties from wild and distant species to broaden the genetic base of cultivated species. For more than 20 years, since Eshed and Zamir identified QTLs related to fruit solids, plant height, and green fruit weight in tomato (Alpert and Tanksley, 1996), CSSLs have received extensive attention. At present, more than 100 sets of CSSLs have been constructed in 20 crops. At the same time, a large number of QTLs were identified by these CSSLs, and fine mapping and gene mining of these QTLs was carried out (Tian et al., 2006;Gu et al., 2015;Wang et al., 2016b;Xin et al., 2016;Zhang et al., 2019).
As high-quality genetic populations, the above CSSLs can help to understand the genetic basis of crop traits and mining important candidate genes. At the same time, marker density has become the limiting factor in gene detection efficiency. Resequencing technology can provide high-density single-nucleotide polymorphism (SNP) markers. Employing SNP data to divide bin markers in the whole population is a feasible marker strategy that can detect QTLs more effectively than traditional methods. At present, it has been employed for fine mapping of QTLs for yield traits in rice, maize, and other crops (Fonceka et al., 2012;Chen et al., 2014;Gu et al., 2015;Wu et al., 2020).
Annual wild soybean is considered to be the wild ancestor of cultivated soybean (Broich and Palmer, 1980). In the process of domestication and improvement, pods, seeds, and other important yield traits have been under constant selection. At present, a large number of seed-pod trait QTLs have been detected. For example,  carried out QTL analysis and locus interaction effect analysis on the number of threeseed pods in a high-generation recombinant inbred line by a SLAF-seq high-density genetic map, and a total of 17 QTLs were detected. Han et al. (2020) performed additive QTL and epistatic QTL mining on the number of seeds, seed weight, number of pods, three-seed pods, and four-seed pods of recombinant inbred lines derived from the crossing of Hefeng 25 and Lmuri 28 by a molecular genetic linkage map, which was constructed by the Bin-Map mapping technique, and 27 related QTLs were obtained. However, only Ln, GA20OX, SWEET10, and PP2C-1 were confirmed to be associated with seed-pod traits (Lu et al., 2016(Lu et al., , 2017Sayama et al., 2017;. Furthermore, the allelic differences between wild and cultivated soybeans for these traits are largely unknown.
In this study, based on the original CSSLs, the genome was resequenced to QTL mapping of seed-pod traits. The purpose of this study was to broaden the genetic base of cultivated soybean, explore the genetic mechanisms of the formation of soybean pod traits, and mine the excellent genes of wild soybean by wild soybean resources. It is expected that the CSSLs will provide excellent breeding and genetic material and support for breeders and geneticists. Revealing the genetic mechanism of traits and mining important candidate genes could provide important theoretical support for the analysis of seed-pod trait formation in soybean.

Field Planting and Management of Soybean
The experimental materials used in this study included 213 CSSLs. The parents were SN14 (recurrent parent) and ZYD0006 (donor parent), constructed with the aid of simple sequence repeats (SSR) markers (Jiang et al., 2020). The materials of CSSLs consist of eight generations (6 BC 3 F 5 , 37 BC 3 F 6 , 61 BC 3 F 7 , 2 BC 4 F 4 , 12 BC 4 F 5 , 22 BC 4 F 6 , 13 BC 5 F 4 , and 60 BC 5 F 5 ). During 2016-2020, the CSSLs were planted at Xiangyang Farm, Harbin, Heilongjiang Province (Harbin, latitude 45°45′N, longitude 126°38′E) in May every year. A randomized complete block design with three replicates was used for 5 years. The lines were grown in a one-row plot with 60 cm row spacing and 6 cm between plants. The rows were 5 m long, with approximately 80 plants per row. The germplasm resource population was planted in Xiangyang Farm, and the planting method is the same as CSSLs. To fine-map the HSW trait, a CSSL (R183) with small HSW was selected as male parent to cross with SN14. Then, the F 1 plants were selfed in Hainan Experimental Station with the false hybrids removed. After selfing, 74 F 2:3 and 144 RHL lines, which comprised the secondary population named R183-F 2:3 and RHL, were harvested. Farmland management of all the above materials follows common agricultural practices. weight (HSW), and seed-pod weight (SPW). Pod-related traits include the proportion of one-seed-pod-numbers (PoOSPN), the proportion of two-seed-pod-numbers (PoTSPN), the proportion of three-seed-pod-numbers (PoTHSPN), and the proportion of four-seed-pod-numbers (PoFSPN). Seed length is the sum of the longest distance parallel to the hilum of 10 seeds. Seed width is the sum of the longest distance perpendicular to the hilum of 10 seeds. The HSW is the weight of 100 seeds of uniform size. Seed-pod weight is the weight of all seeds and pod shells on the whole plant .

Material Resequencing and Bin Marker Acquisition
Genomic DNA was extracted from young leaves of Suinong 14 (SN14), wild soybean ZYD00006, and CSSLs by a modified CTAB method (Wu et al., 2020). Then, the extracted genomic DNA was randomly interrupted by ultrasonic fragmentation, and the sequencing library was constructed by terminal repair, adding an a at 3′ end, adding a sequencing connector, purification, and PCR amplification. After passing a quality inspection, the library was sequenced with an Illumina HiSeqTM sequencing platform, and the sequence read length was 150 bp. By BWA and GATK software, a sequence was compared with the reference genomic Glycine max Wm82.a2.v1 to identify the original SNP (Li and Durbin, 2009;McKenna et al., 2010;Schmutz et al., 2010), and high confidence SNPs were obtained according to the method of Guo et al. (2015).
A developed approach of sliding-window (Huang et al., 2009) was used for CSSLs genotyping with a group of 17 consecutive SNPs and a step size of 1. The position of a SNP that switched one genotype to another consecutive genotype was recorded as the recombination breakpoint. A bin map was constructed according to the recombination breakpoint information. All of the breakpoints along the same chromosome were sorted by physical position.

Method of Date Analysis
Phenotypic data of CSSLs and secondary segregation population (F 2 , F 2:3 , and RHL) have been statistically analyzed by SPSS 26. For example, average value, maximum value, minimum value, etc. The best linear unbiased estimate (BLUE) can more effectively eliminate the deviation caused by the environment and improve the accuracy of QTL mapping. The BLUE was calculated by BLUE calculation option in the AOV module of the ICIMapping 4.2 software. The parameters used were each line of CSSLs as fixed factor, years and repetition as random factors.

QTL Mapping of Pod-Seed Traits
Bin markers were used for QTL mapping of the seed-pod-related traits in CSSLs, and a likelihood ratio test based on stepwise regression model (RSTEP-LRT-ADD) with ICIMapping 4.2 software was used to determine QTLs. Parameters setting was "By condition number" = −1,000 (equivalent to deleting duplicate markers), "PIN" = 0.001 (PIN: the largest-value for entering variables in stepwise regression of residual phenotype on marker variables), and LOD value ≥2.5 (Wang et al., , 2016aLi et al., 2007Li et al., , 2008. The linear model of the RSTEP-LRT method was as follows: where i = 0 (for the recurrent parent), 1, 2, …, n, b 0 is the intercept, b j (j = 1, …, t) is the partial regression coefficient of phenotype on the jth marker, which represents QTL additive effect on each segment, x ij is the indicator variable for the jth marker in the ith CSSL, which is equal to −1 if the marker type is the same as in P 1 (the recurrent parent) and 1 if the marker type is the same as in P 2 (the donor parent), and e i is the random experimental error following a normal distribution. DNA was extracted from the secondary segregated population by the CTAB method. SSR markers were used as molecular markers (the marker obtained from Soybase). 1 After polymerase chain reaction (PCR) of Panaud, the genotype of secondary segregated population was then detected by 8% PAGE separation (Doyle and Doyle, 1990). Genetic map is constructed by MAP module in ICIMapping 4.2 software. QTL mapping employed the ICIM-ADD model with 1,000 permutations calculation of ICIMapping 4.2 software. And the threshold of ICIMapping for QTL detection was chose p ≤ 0.05. The epistasis of the seed-pod-related traits by the RSTEP-LRT-EPI model with PIN = 0.000001calculation of ICIMapping 4.2 software . And the threshold of ICIMapping for epistasis QTL detection were set chose LOD > 10.

qRT-PCR Analysis of Candidate Genes
The seed of qRT-PCR material was sampled according to different developmental stages of soybean seeds. The Cotyledon Stage (COT, 10 days after flowering (DAF)) and Early Maturity Stage 1(EM1, 20 DAF) mainly carry out cotyledon development and endosperm assimilation. Then, the seed size and weight are increasing at Early Maturity Stage 2(EM2, 30 DAF) and Mid Seed Maturity (MM, 50 DAF). 2 During the COT, EM1, EM2, and MM of seed in parents SN14, CSSL-R183, and wild soybean ZYD00006 (Severin et al., 2010;Collakova et al., 2013;Lu et al., 2016), soybean seed RNA was extracted by TRIzol (Invitrogen, 15596-026, Carlsbad, CA, United States). qRT-PCR was carried out on a Roche LightCycler 480 II fast real-time PCR system using SYBR qPCR Mix (Vazyme, Q711, Vazyme Biotech, Nanjing, China). Actin (GLYM18G52780) was used as the internal reference (Gu et al., 2017), and the relative amount of gene expression (target gene/actin) was calculated by the 2 ΔCT method. Premier 5.0 3 was used to design qRT-PCR primer sequences for candidate genes.

Haplotype Analysis
A total of 129 germplasm resources from China were used for haplotype analysis due to their rich genetic variation and large 1 https://soybase.org/ 2 https://soybase.org/amigo/amigo/term/SOY:0001283#display-lineage-tab 3 http://www.premierbiosoft.com/primerdesign/index.html range of HSWs. The genomic sequences of the candidate genes were proposed on the Phytozome13 website, 4 including the 5′UTR, upstream 2000 bp promoter region sequence, CDS, intron sequence, and 3′UTR sequence information. The sequence information for the candidate gene and the genomic sequencing information from the resequencing of 129 soybean germplasm resources were subjected to local BLAST analysis to obtain SNP information from the candidate genes in the germplasm resource population (unpublished). In this study, Dnasp5.0 software was used to analyze the haplotype distribution of candidate gene SNP sequences in germplasm resource populations, and excellent haplotypes (the number of cultivars belonging to the haplotype exceeding 5% of the total) were screened out. 5 Using the Haps Format module in Haploview 4.2 software (Cambridge, Massachusetts, United States) to analyze the excellent haplotype sequence information from candidate genes, ANOVA of candidate gene haplotypes and phenotypes was performed by SPSS 26.0 software (Armonk, New York, USA), to determine the effect of each haplotype on phenotype. Plant cis-acting regulatory DNA elements (PLACE) 6 and plant promoters and cis-acting regulatory elements (Plant Care) 7 were used to query the function of plant promoter elements.

Resequencing of CSSLs (Parents) and Division of Bin Markers
The sizes of the sequenced data from SN14, ZYD00006, and the 213 CSSLs were 30.08 Gb, 33.38 Gb, and 1410 (Schmutz et al., 2010), and the results showed that the average alignment efficiency of SN14, ZYD00006, and the CSSLs was more than 90%. The average genome coverage depth and percentage of parents were more than 27× and 95% (covering at least 1×), respectively. The average genome depth and percentage of the CSSLs were more than 5.79× and 95.74% (covering at least 1×), respectively (Tables 1, 2). At the same time, 580,524 SNPs were identified between SN14 and ZYD00006, which indicated that there were genetic differences between cultivated soybean SN14 and wild soybean ZYD00006. These SNPs in the parents of the CSSLs also facilitate further fine mapping of QTLs.
The 580,524 SNP was divided into 3,780 bin markers. The range of bin marker length was 20 Kb-22.91 Mb, and the average marker length was 225.12 Kb. 79.26% of bin markers were less than 0.2 Mb in length, and 111 markers were more than 1 Mb in length ( Figure 1A). There were 6 bin markers greater than 10 Mb, which were distributed on Chr. 04 (Bin2350), Chr. 08 (Bin4681), Chr. 09 (Bin4865), Chr. 12 (Bin6269), Chr. 18 (Bin9389), and Chr. 19 (Bin9688). At the same time, there were differences in the number and density of bin markers on 20 soybean chromosomes ( Table 3). The highest and lowest bin marker densities were on chromosome 9 (2.33 markers/Mb) and chromosome 3 (5.50 markers/Mb), respectively. The average marker density was 4.04 markers/ Mb, and the average distance between bin markers was 25.45 Kb. These CSSLs carried 2441 introgressed segments (including 459 introgressed segments with a single bin marker). The length of these introgressed segments ranged from 0.02 Mb to 51.74 Mb. Among the above introgressed segments, 60.43% of the introgressed segments were less than 1 Mb, and large segments greater than 10 Mb accounted for 7.87% ( Figure 1B). The CSSLs many large introgressed fragments, so they need to be further expanded.

Phenotypic Variation and Statistical Analysis of CSSLs
The phenotypic data of CSSLs from 2016 to 2020 were statistically analyzed and the BLUE was calculated.  Table 2). The above results showed that the phenotypic variation in the CSSLs was plentiful, indicating that the CSSLs were suitable for genetic mapping research. The correlation between seed-pod-related traits was analyzed. PoOSPN and PoTSPN were positively correlated;  PoTHSPN and PoFSPN were positively correlated; PoOSPN were also negatively correlated with PoTHSPN and PoFSPN; and PoTSPN were also negatively correlated with PoTHSPN and PoFSPN. SL was positively correlated with SW, SL and SW were positively correlated with HSW. SL, SW, and HSW were significantly positively correlated with SPW. The above results showed that the pod-seed-related traits of soybean were not isolated, and the traits were closely related (Supplementary Table 3).

Screening Materials for Further Construction of CSSLs
Based on the information of sequencing and introduction of donor fragments, 129 lines which can be achieve the maximum coverage of CSSLs were selected from 213 CSSLs as the basis for future construction of CSSLs (Figure 3). These 129 CSSL carried 1772 introgressed segments (including 359 introgressed segments with a single bin marker; Figure 4). The average number of introgressed segments per chromosome in these lines was 88.6. The length of these introgressed segments ranged from 0.02 Mb to 47.14 Mb, and their average length was 3.3 Mb. Among the above introgressed segments, 62.19% of the introgressed segments were less than 1 Mb, 24.49% were between 1 and 5 Mb, and large segments greater than 10 Mb accounted for only 7.45%. The total length of introgressed segments in the CSSLs was 5589.34 Mb, which was 5.71 times larger than the soybean genome and covered 89.65% of the genome of wild soybean (Table 4). In the future, we will further expand the CSSLs on this basis of 129 lines to complete the construction of whole-genome substitution lines.    Table 5). The PoOSPN and SPW were no QTL detected steadily for more than 2 years. Compared with total QTLs mapped for each trait, 23, 7, 1, and 1 QTLs were repeatedly detected in 2, 3, 4, and 5 traits, respectively. Among the QTLs that were detected steadily for more than 2 years, QTLs with the highest proportions of phenotypic variation explained were qPoTSPN13-1, qPoTHSPN13-1, qPoFSPN20-3, qSL17-1, qSW17-1, and qHSW17-3, respectively ( Table 5). qPoTSPN13-1 and qPoTHSPN13-1 are the same locus detected on chromosome 13, and the proportions of phenotypic variation explained were 20.48 and 16.14%, respectively. qPoFSPN20-3 is a QTL detected on chromosome 20, and the proportion of phenotypic variation explained was 30.69%. qSL17-1 and qSW17-1 are the same locus detected on chromosome 17, and the proportions of phenotypic variation explained were 28.76 and 16.54%, respectively. qHSW17-3 is a QTL detected on chromosome 17, and the proportion of phenotypic variation explained was 37.58%. When comparing QTLs related to different traits that were stably detected for more than 2 years, we found that two QTLs, qPoTSPN18-1  and qSL18-1, were detected at locus Bin9163; three QTLs, qPoTSPN10-1, qPoFSPN10-1, and qSW10-2, were detected at locus Bin5399; two QTLs, qPoTSPN13-1 and qPoTHSPN13-1, were detected at locus Bin6769; three QTLs, qPoTSPN20-3, qPoTHSPN20-1, and qPoFSPN20-3, were detected at locus Bin10791; and the locus Bin8866 was detected as qSL17-1, qSW17-1 and qHSW17-4 ( Figure 5). Compared with the QTL detected only once, the QTL repeatedly mapped between different traits or between different years of the same trait has higher feasibility. The above QTLs will be the focus of further research.

Establishment of a Secondary Segregating Population and Fine Mapping
R183 was a material in CSSLs with extremely small HSW; several HSW QTLs stably detected more than 2 years were located in the introgressed segment of R183 (bin2808, bin5553, and bin6188). The secondary segregating population was constructed by backcrossing R183 (12.74 g) with SN14 (19.73 g; Supplementary Figures 2A,B) to identify candidate genes regulating soybean HSW. Forty-five pairs of SSR markers with polymorphisms between parents were screened. Preliminary mapping of QTLs employing the secondary segregating population R183-F 2 and the 45 SSR markers detected 5 QTLs (Supplementary Figure 2C; Table 6). Among them, qHSW-12-1 was located between two flanking markers, M371 (7.27 Mb) and M570 (11.68 Mb). qHSW-12-1 had the highest proportion of phenotypic variation explained (21.42%). Next, fine mapping of qHSW-12-1 was performed ( Figure 6A).
One heterozygous plant in the target interval was screened from the secondary segregating population R183-F 2 , and 74 R183-F 2:3 plants were obtained to narrow the target interval. Fine mapping of qHSW-12-1 was performed by the secondary segregating population R183-F 2:3 and 5 SSR markers in the target interval. qHSW-12-1 was mapped between two flanking markers, M371 (7.27 Mb) and M454 (8.50 Mb; Supplementary Figures 2D,E). The fine mapping QTL interval was larger, so individuals with heterozygous target segments were selected to construct the RHL population, and a total of 144 individuals were obtained. Ten SSR markers were employed to increase the density of primers in the target interval. Finally, qHSW-12-1 was narrowed down to a 64 kb region between M371 (7.27 Mb) and M374 (7.33 Mb; Figure 6B; Supplementary Figures 2F,G).
According to the Williams 82 soybean reference genome sequence in the Phytozome database (See Footnote 4), the genome interval of qHSW-12-1 contains a total of 5 genes (Glyma.12G088700, Glyma.12G088800, Glyma.12G088900, Glyma.12G089000, and Glyma.12G089100; Figure 6B). Glyma.12G088700, a transcription factor related to drought resistance and salt tolerance (Liao et al., 2008;, showed no difference between the parents in the sequences of the upstream 2000 bp promoter or the CDS region. Glyma.12G088800 has no genetic annotation, and quantitative real-time PCR (qRT-PCR) showed that the gene was not expressed during seed development periods. Glyma.12G089000 (annotated as transmembrane 9 superfamily member 5) is a gene that has not been studied at present. qRT-PCR showed that there was a significant difference in gene expression between parents at the COT and EM1 stages, but there was no significant difference at the EM2 and MM stages. Glyma.12G089100 is also called GmPHR20. Its homologous gene (GmPHR25) in soybean is an important regulator of the phosphorus signal network (Xue et al., 2017), and its homologous gene (AT3G04030) in Arabidopsis may be involved in nitrogen reuse (Nakano et al., 2017). However, qRT-PCR showed that there was no difference in the relative expression of the gene between parents.  There were 7 SNP differences in the promoter region and 2 SNP differences in the CDS region in Glyma.12G088900 between parents (Supplementary Figures 3-5). qRT-PCR showed that there were extremely significant differences in the gene during the seed development periods of both parents (Figure 7A; Supplementary Table 5). Dnasp5.0 software was used to analyze the haplotype distribution of Glyma.12G088900 SNP sequences in germplasm resource populations. There was an extremely significant difference in HSW between two excellent haplotypes, Hap1, with the same genotype as R183, and Hap2, FIGURE 5 | The QTLs were detected steadily for more than 2 years. The QTLs detected steadily for 2 years were yellow; The QTLs detected steadily for 2 years were green; The QTLs detected steadily for 5 years were red.

Analysis of Epistasis Between Loci
The epistatic effects were analyzed by the BLUE values of 8 traits in 5 years. The results showed that 955 pairs of epistatic QTLs were obtained from the other 6 traits except PoFSPN and SPW (Figure 8; Supplementary Table 7). Among these epistatic QTLs, the LOD values were between 10.01 and 36.08, and the proportions of phenotypic variation explained were between 0.11 and 12.85%. There were 171, 79, 6, 38, 282, and 379 pairs of epistasis in HSW, PoOSPN, SL, SW, PoTSPN, and PoTHSPN, respectively. Among them, the maximum proportions of phenotypic variation explained were 0.63, 0.46, 12.85, 1.60, 0.53, and 0.35%, respectively, and the minimum proportions were 0.14, 0.28, 5.52, 1.12, 0.16, and 0.11%, respectively. Except for SL, the proportions of phenotypic variation explained of epistasis about other traits were lower, indicating that the main effect QTL played an important role in the regulation of target traits. In HSW, 19 pairs of epistasis occurred between two major QTL, the other occurred between a major QTL and a locus, and a total of 129 pairs played a positive role in increasing HSW. qHSW17-3 not only had epistatic interaction with 5 major QTL (qHSW3-2, qHSW5-3, qHSW4-1, qHSW18-2, qHSW20-3), but also had epistatic interaction with 110 motifs, indicating that the QTL may play an important role in the regulation of HSW.
Seventy-eight pairs of epistasis related to PoOSPN occurred between one major QTL and one locus, only one pair of epistasis occurred between two major QTL. Among them, 59 pairs played a positive role in reducing PoOSPN. All the above 79 pairs of epistasis were related to qPoOSPN14-1, indicating that qPoOSPN14-1 played an important role in regulating PoOSPN.
Six pairs of epistasis related to SL occurred between one major QTL and one locus, in which the pairs of epistasis between qSL16-1 and 2 loci had negative effects on SL increase and the pairs of epistasis between qSL17-1 and 4 loci played positive roles in SL increase. The results showed that qSL17-1 may play an important role in regulating SL.
Thirty-eight pairs of epistatic effects related to SW occurred between a major QTL and a locus, in which 20 pairs of epistasis had negative effects on the increase of SW. There were epistatic interactions between qSW1-2 and 34 loci, indicating that qSL17-1 may play an important role in regulating SW.
Eighty-eight pairs of epistasis related to PoTSPN occurred between two major QTL, and the other occurred between one A B FIGURE 6 | QTL fine mapping. (A) The qHSW-12-3 locus was detected on chromosome 12. Positional mapping narrowed the qHSW-12-3 locus to a 64.4 kb region between M371 and M374. The numbers above the bars indicated the number of recombinants. (B) Fine-map the genotype of the population and the significance of the corresponding phenotype. Gray indicates Suinong14 genotype, black indicates ZYD00006 genotype. a, b: Different letters represent significant differences between each other at the 0.01 level.
major QTL and one locus. A total of 261 pairs of epistatic played positive roles in reducing PoTSPN. qPoTSPN11-1, qPoTSPN13-1 and qPoTSPN20-3 not only had epistatic interactions with 6, 4, and 4 major QTLs, respectively, but also had epistatic interactions with 14, 27, and 28 loci, respectively, indicating that the three QTLs may play important roles in regulating PoTSPN.
One hundred sixty pairs of epistasis related to PoTHSPN occurred between two major QTLs, and the other occurred between one major QTL and one locus. qPoTHSPN11-1, qPoTHSPN13-1, and qPoTHSPN20-1 had epistatic interactions not only with 4, 2, and 4 major QTLs, respectively, but also with 51, 56, and 45 loci, respectively, indicating that the three QTLs may play important roles in regulating PoTHSPN.

DISCUSSION
Wild soybean was domesticated approximately 6,000-9,000 years ago (Zhou et al., 2015), and is recognized as the wild ancestor of cultivated soybean (Zhang et al., 2018). In the process of domestication, due to purposeful selection by human beings, the genetic diversity of cultivated soybean became much lower than that of wild soybean (Hyten et al., 2006;Lam et al., 2010). In this study, CSSLs were constructed between cultivated soybean SN14 and wild soybean ZYD00006. Through multigeneration backcrossing, the introduction of the wild soybean genome greatly enriched the genetic base of cultivated soybean, which led to rich phenotypic variation in the CSSLs. For example, the maximum HSW is 25.69 g, the minimum is 9.83 g; the maximum SPW is 38.68 g, the minimum is 11.62 g; and the maximum PoTSPN is 59.06%, the minimum is 15.08% (Figure 2; Supplementary Table 1). In this study, we also employed SNP markers to divide the soybean genome into 3,780 bin markers. High-density bin markers can improve the efficiency and accuracy of QTL mapping, which makes it possible to perform fine mapping of QTLs and gene mining.
In the next step of research, 129 lines covering 89.65% of wild soybean genome were selected from CSSLs and hybridized with SN14. MAS was used to shorten the introgressed fragments, purify the genetic background of the CSSLs and improve the genome coverage of wild soybean. Finally, we will complete the construction of whole-genome substitution lines so that they have greater utilization value.
In order to increase the accuracy of QTL mapping, we integrated the phenotypes in different years to calculate the BLUE value of each trait. BLUE, derived from the linear mixed model, has been popularly used to estimate animal and plant breeding values for a few decades (Kim, 2020). BLUE can more effectively eliminate the deviation caused by the environment and improve the accuracy FIGURE 8 | Epistasis analysis of CSSLs. Epistasis of seed-pod traits that have repeatedly appeared over the years. Dark blue represents the epistatic effect among HSW loci, red represents the epistatic effect among PoOSPN loci, black represents the epistatic effect among SL loci, yellow represents the epistatic effect among SW loci, green represents the epistatic effect among PoTHSPN loci, gray represents the epistatic effect among PoTSPN loci. of QTL mapping. Muqaddasi et al. (2017) studied anther extrusion of wheat accessions by genome-wide association studies with BLUE of anther extrusion. Then, they predicted the candidate genes which could be directly or indirectly involved in anther extrusion.  employed the BLUE of soybean plant height and the number of nodes on the main stem for QTL mapping, and predicted two candidate genes for simultaneous regulation of these two traits. In addition, we selected appropriate method to improve the accuracy of QTL mapping. Because most CSSLs used for QTL mapping are non-idealized CSSLs (There is only one import fragment in each CSSL.; Li et al., 2015;Bux et al., 2019). There are redundant markers, and there may be multicollinearity between markers. The goal of QTL mapping is to explore the impact of different markers on target traits and quantitatively describe their impact on target traits. Multicollinearity between markers will reduce the positioning ability and calculation accuracy, so multicollinearity must be controlled. Wang et al. (2006) developed RSTEP-LRT model for QTL mapping of non-idealized CSSLs. Firstly, the model will evaluate the multicollinearity between different markers and delete redundant markers. However, there is no theoretical method to find an appropriate threshold to calculate whether the multicollinearity is high. In this study, we used the condition number −1,000 as the threshold of multicollinearity, which was one preset value in the software. It is equivalent to deleting duplicate markers, to reduce the adverse effect of multicollinearity on QTL mapping.
QTL mapping is based on the difference in natural alleles among parent varieties (Ma et al., 2019), and its mapping accuracy is limited by the mapping population. CSSLs can eliminate background interference during QTL mapping and maintain high mapping accuracy when QTLs are located on very small fragments. In recent years, CSSLs have been widely used in QTL mapping and candidate gene mining (Balakrishnan et al., 2019). Huo et al. (2017) used CSSLs to fine map the number of grains per panicle and verify the function of the candidate gene NOG1. Yang et al. (2017) fine mapped soybean HSW, protein content, and oil content using CSSLs. The bidirectional sucrose transporter gene Glyma.15 g049200 was predicted to be a pleiotropic regulatory gene that regulates HSW, protein content, and oil content. In our study, QTL mapping of seed-pod traits in wild soybean was carried out by CSSLs of wild soybean. Eight QTLs related to HSW were identified for more than 2 years, of which four QTLs coincided with previous mapping results (Hyten et al., 2004;Gai et al., 2007;Kuroda et al., 2013), namely qHSW9-2, qHSW11-2, qHSW17-3, and qHSW17-3. Five-seed length QTLs and four seed width QTLs were reported steadily for more than 2 years, and two of them were seed width QTLs located on chromosome 1 for the first time. qPoTSPN20-3, qPoTHSPN20-1, and qPoFSPN20-3 are the same locus detected on chromosome 20. This is the only QTL detected that has been stable for more than 2 years and is related to the POTSPN, PoTHSPN, and PoFSPN. Compared with the results of previous studies, we found that the QTL detected on chromosome 20 coincided with the cloned Ln (Sayama et al., 2017), and Ln had been demonstrated to have a polymorphism to control the leaf shape and pod number for soybean. The above results reflect the accuracy and effectiveness of this set of CSSL materials in determining the positions of QTLs.
To explore the excellent alleles regulating HSW, we constructed a secondary segregation population by backcrossing CSSL-R183 (which had an extremely small HSW) with SN14. The HSW QTL was mapped to the 64.4 Kb region of chromosome 12 by two fine mapping. Glyma.12G088900 was preliminarily identified as a candidate gene by sequence analysis of genes in the 64.4 kb region, relative expressivity analysis at important developmental stages, and haplotype analysis. Glyma.12G088900 belongs to the F-box/LRR family of proteins, which can regulate a variety of activities, such as regulating plant floral organ development, responding to biotic and abiotic stresses, and regulating growth, development and crop yield (Chen et al., 2013;Jia et al., 2015;Baute et al., 2017;He et al., 2017). At the same time, the homologous gene in Arabidopsis (AT1G67190) can interact with SKP1-like protein to form the SCF complex and regulate grain development, in which ASK1 and ASK2 have the function of regulating embryogenesis and seed development (Liu et al., 2004;Kuroda et al., 2012).
Many studies have shown that epistatic effects were important factors affecting agronomic traits of soybean and important genetic components of quantitative trait variation, and the epistatic effects between major QTLs had very important effects on quantitative traits (Jiang et al., 2018). Qi et al. (2017) used RIL and CSSLs to carry out QTL mapping and epistatic analysis of soybean oil traits, and a total of 11 QTLs and 4 pairs of epistatic QTLs were detected, which provided important information for soybean molecular-assisted breeding. Balakrishnan et al. (2019) used CSSLs to carry out QTL mapping and epistatic analysis. A total of 22 QTLs and 17 pairs of significant epistatic QTLs for 11 traits were detected, which provided valuable resources for gene discovery and yield improvement. In this study, 955 pairs of epistatic QTLs were detected by CSSLs to analyze seed-pod traits of soybean. Among them, qHSW17-3, qPoOSPN14-1, qPoTSPN11-1, qPoTSPN13-1, qPoTSPN20-3, qPoTHSPN11-1, and qPoTHSPN20-1 have epistatic interactions not only with other major QTLs, but also with a large number of loci. But the proportions of phenotypic variation explained of above epistasis were low, the highest is only 0.62%. However, when a large number of epistatic effects were aggregated, it was still an important part that cannot be ignored in actual production. For breeding programs to improve seed-pod traits, in addition to the most important process of aggregating reliable QTLs (especially the main QTLs), attention should also be paid to the loci and interaction effects that interact with main QTLs, then the QTLs can make the greatest contribution to target traits (Jiang et al., 2011). In addition, bin5858 (qPoTSPN11-1, qPoTHSPN11-1) had epistasis interactions with bin6769 (qPoTSPN13-1, qPoTHSPN13-1) and bin10791 (qPoTSPN20-3, qPoTHSPN20-1) in PoTSPN, PoTHSPN, respectively. Interestingly, the two pairs epistasis resulted in a decrease in PoTSPN, while an increase in PoTHSPN. The results showed that seed-pod traits were often closely related, and the same locus can not only directly affect different pod characters, but also indirectly regulate them through epistasis (Fang et al., 2017;Yang et al., 2021). The loci controlling several pod traits and the loci with a large number of epistasis found in this study can provide valuable information for molecular-assisted breeding of soybean.

CONCLUSION
In this study, 580,542 SNPs were identified between SN14 and ZYD00006, and were divided into 3,780 bin markers. A total of 3,780 bin markers were used for QTL mapping of seed-podrelated traits in CSSLs. A total of 170 QTLs were detected, of which 32 QTLs appeared stably in different years. A secondary segregation population was constructed for HSW, and qHSW-12-1 was mapped to the 64.4 Kb region of chromosome 12. qRT-PCR showed that there were significant differences in the gene during seed development in the two parents. Further analysis showed that there was a difference in Glyma.12G088900 in the promoter region and an amino acid difference in the coding region, which may affect HSW. A total of 955 pairs of epistatic sites were detected in epistatic analysis. Among them, 161 QTL mapped by CSSLs have epistatic effect. The construction of wild soybean CSSLs provides excellent genetic materials for broadening the genetic base of cultivated soybean and mining important candidate genes, and the genetic study of related traits lays a theoretical foundation for the study of seed-pod traits in soybean.

DATA AVAILABILITY STATEMENT
The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found in the article/Supplementary Material.

AUTHOR CONTRIBUTIONS
CL, HJ, and QC conceived the study and designed and managed the experiments. MY, CL, and ZQ provided plant lines. HZ, FC, and RW performed trials and collected data. ZZ, RZ, DX, and HZ completed statistical analyses of phenotypic data. HZ, LH, and JX contributed to writing the paper. All authors contributed to the article and approved the submitted version.


s2fieldsofstudy:
Agricultural And Food Sciences
Medicine

text:
Cow's milk and egg protein threshold dose distributions in children tolerant to beef, baked milk, and baked egg

Abstract Background The use of eliciting doses (EDs) for food allergens is necessary to inform individual dietary advice and food allergen risk‐management. The Eliciting Dose 01 (ED01) for milk and egg, calculated from populations of allergic subjects undergoing oral food challenges (OFCs), are 0.2 mg total protein. The respective Eliciting Dose 05 (ED05) is 2.4 mg for milk and 2.3 mg for egg. As about 70% children allergic to such foods may tolerate them when baked, we sought to verify the EDs of that subpopulation of milk and egg‐allergic children. Methods We retrospectively assessed consecutive OFC for fresh milk and egg between January 2018 and December 2020 in a population of baked food‐tolerant children. Results Among 288 children (median age 56 ‐ IQR 36–92.5 months, 67.1% male) included, 87 (30.2%) returned positive OFC results, 38 with milk and 49 with egg. The most conservative ED01 was 0.3 mg total protein (IQR 0.03–2.9) for milk and 14.4 mg total protein (IQR 3.6–56.9) for egg. The respective ED05 was 4.2 (IQR 0.9–19.6) mg for milk and 87.7 (IQR 43–179) mg for egg. Such thresholds are, respectively, 1.5 (milk ED01), 1.75 (milk ED05), 72 (egg ED01), and 38.35 (egg ED05) times higher than the currently used thresholds. Conclusions The subpopulation of children allergic to milk and egg, but tolerant to baked proteins, displays higher reactivity thresholds than the general population of children allergic to milk and egg. Their risk stratification, in both individual and population terms, should consider this difference. In baked milk‐tolerant children, milk causes reactions at lower doses than egg in our group of egg‐tolerant children. This could be associated with the relative harmlessness of egg compared with milk in the determinism of fatal anaphylactic reactions in children.

| INTRODUC TI ON
Allergic reactions to milk and egg proteins are the most prevalent food allergies in European children. 1,2 They are a common cause of severe life-threatening reactions in pediatric allergic patients and could persist throughout life, placing a heavy burden on the quality of life, the development of other allergic diseases, and growth. [3][4][5] In the most sensitive individuals, trace amounts of milk and egg allergens may provoke allergic reactions, even from incidental contacts as opening packages, inhaling vapors from cooking, use of shared utensils, and kissing the lips of someone who has eaten the offending food. 6,7 Such reactions may be life-threatening, in particular for cow's milk. This food triggers more than 25% cases of fatal anaphylaxis in children, while egg has only exceptionally been associated with mortality for food allergy. 8 The primary and safest strategy for managing milk/egg allergy is strict avoidance of the causal food. To implement it, clinicians personalize their advices to patients informing the suggestions to their threshold of reactivity at oral food challenges (OFCs) among other clinical considerations. The population of milk-and eggallergic children includes a spectrum of more and less sensitive individuals. The former may also react to milk and egg in baked products as muffins or biscuits; less sensitive children tolerate baked foods. 9 In addition, the 10-20 percent milk-allergic children reacting to beef are considered having a more severe and persistent Cow's Milk Allergy (CMA). 10,11 When OFC-negative to baked milk and egg, children are not requested to avoid such foods. 12 Foodallergic children with a high reactivity threshold may not be requested to avoid precautionary-labelled prepackaged foods. 13,14 The use of precautionary allergen labelling (PAL) is based on population-based threshold dose distributions for the different food allergens. 15,16 To this end, data from cohorts of subjects assessed in oral food challenges (OFCs) are used to calculate eliciting doses (EDs) and to specify the relation between EDs and symptoms. 17,18 Based on such thresholds, the Voluntary Incidental Trace Allergen Labelling (VITAL), a scheme developed by the Australian and New Zealand food industries, established reference doses of 0.2 mg total protein (Eliciting Dose 01) for milk and egg, 16,19 and an Eliciting Dose 05 between 0.5 and 2.4 mg protein for both foods. 16,20,21 Such thresholds have been calculated on populations of patients allergic to both raw and baked foods. As raw and baked food are not likely equivalent, we aimed to explore the EDs of the subset of children allergic to milk and egg, who tolerate beef and baked proteins. We posed the hypothesis that children allergic to native products, but not to baked foods, may have a higher reactivity threshold that may exempt them from compliance with the current PAL.

| Study design
At the Bambino Gesù Allergy Department, patients with suspected milk/egg allergy are periodically exposed to oral food challenges (OFCs) until they get tolerance. For safety reasons, patients are excluded from OFC if an anaphylactic reaction had occurred in the last 6 months for children of 0.5-5 years, 12 months for children of 6-12 years, and 2 years for patients aged over 13 years. 22 Also not In this setting, we retrospectively assessed consecutive children who underwent OFC for fresh milk and egg between January 2018 and December 2020. For those who underwent repeated OFCs to the same food, we included only the first procedure for each type of food.

G R A P H I C A L A B S T R A C T
The clinical suspicion was placed based on the clinical history, and the entry symptoms were classified on parents' descriptions reported in the personal medical records. We only considered children with a personal history of acute IgE-mediated allergic reactions (developing within 1-2 h after food intake) from mild-to-severe systemic reactions (anaphylaxis). Before the OFC, parents were exhaustively informed about the risks of the procedure and gave their written consent. On the day, we administered the OFC, the sensitivity to the respective food was tested at skin prick test (SPT) and spe-

| Symptom grading and threshold data
Our trained staff performed OFCs in hospital with close supervision and immediate availability of emergency treatment, according to guidelines. The same clinicians observe and value the oral food challenges' result, and the doses were administered until objective symptoms appeared. 25,26 Briefly, the reactive symptoms classified patients in five groups, from subjective symptoms (nausea, abdominal pain, pruritus, oral allergy syndrome) to systemic reactions ( Table 2). 27 For this study, we designated the initial point at which objective symptoms occur at a specific dose as the Lowest Observed Adverse Effect Level (LOAEL) and the highest dose that does not lead to objective symptoms as the No Observed Adverse Effect Level (NOAEL).
According to previous studies, we choose an Interval-Censoring Survival Analysis (ICSA) approach to our data. This method is necessary when the exact dose that provokes a reaction in a patient is unknown, but it falls into a particular interval. Specifically, as previously described, 28 if a patient had an objective reaction at the first dose, left censoring occurred, and the NOAEL was set to zero with the LOAEL set as that first dose. If a patient did not experience an objective reaction after the 7th dose, the cumulative given dose was designated as NOAEL. In this case, we considered the data to be right-censored, and the LOAEL set to infinity. If the patient was unable or unwilling to continue the challenge because of the subjective symptoms, we discontinued the challenge and considered the final administered dose right censored. In all other cases, intervalcensoring occurs bounded by the NOAEL and LOAEL.

| Statistical analysis
We correlated the individual threshold to age and sex of patients, history of clinical manifestations, sIgE level results, or SPTs wheal size using multiple logistic regression (SPSS package, BM Corporation). A p-value of <.05 was considered statistically significant; sIgE values >100 kU/L received a designated value of 101 kU/L.
We also used SAS 15.1 (SAS Institute) and its LIFEREG procedure to fit parametric models to the interval-censored data. We considered Log-normal, log-logistic, and Weibull models to fit these data and to extrapolate EDs (ED01, ED05, ED10, ED25, and ED50), which are the doses predicted to elicit allergic reactions in 1%, 5%, 10%, 25%, and 50% of milk and egg-allergic patient, respectively. 29 The goodness-of-fit was assessed by the Anderson-Darling test.
Confidence intervals (CI) were calculated and added.
Step   Table S1 in this article's Online Repository).

| Eliciting doses
We fitted the data with log-normal, log-logistic, and Weibull distributions models (Figures 2A,B)

| DISCUSS ION
This is the first study aimed to evaluate the thresholds of children allergic to native egg and milk only. Our patients were tolerant to beef, baked milk, and baked egg. This population, which can be identified by specific challenges to baked foods, 30 represents the majority of milk/egg allergy sufferers. As they display a shorter duration of the condition, and present with less complex clinical pictures, they are considered those with the less severe forms of milk and egg allergy. 31 The prevalence of positive OFC in our center was 26.1% and 34.5% for milk and egg, respectively. While a similar prevalence was reported in children with OFCs positive to egg (33%) and milk (32%), 32 other caseloads showed higher prevalence rates, both for milk and egg. [33][34][35][36][37] Our threshold values are similar to those found in VITAL for milk  Table 2) VITAL, consistently with the fact that they are allergic to native egg only and that no statistical correlation between patients' individual threshold and their age, gender, sIgE, and/or SPT levels was found.
Our milk allergic children also showed a higher reactivity threshold than that of VITAL both at level 01 and at level 05, but differences were much smaller. Again, this is unlikely to be due to different food doses at challenge. Studies with milk OFCs in general use inferior cumulative amounts of milk, with starting doses similar [48][49][50] or higher, 51,52,53 to our 3 mg milk protein.
These findings echo those of other studies for milk, but not for egg. For instance, in a single-center pediatric experience using a cumulative dose amount of protein of 2190 mg, corresponding to about 63 ml of milk and 1/3 egg, the ED05 were set at 1.07 for milk and at 1.51 for egg. 54 The only study comparing thresholds of baked-tolerant vs. baked-sensitive egg-and milk-allergic patients found ED10 and ED50 values similar to ours for milk, but much lower for egg 55 (Table 4). We recently published tolerance threshold data of children with severe baked milk and baked egg allergy, using a different OFC scheme than the one used in this study. The mean thresholds of reactivity we found in that context were 116.3 (±107.6) and 128.3 (±96.7) mg protein for milk and egg, respectively; the EDs were not calculated. 22 A further study found much lower levels ED05 and ED10 for milk. 56 These differences could be due to different selected populations, difference in mean age, and difference in the method of administration of the OFC. In any case, should our results be confirmed by other studies they could explain why, in spite of a similar prevalence in the population, milk continues to be more represented than egg in the caseloads of mortality due to food anaphylaxis in children. 8,59 It is licit to infer that lower thresholds may induce parallel higher probability of severe reactions. Despite their low ED01, we found a low frequency In most studies, EDs for cow's milk and egg proteins have been associated to the age of patients, the severity of allergic reactions or the OFC protocol used, such as the number of doses and the total amount of food allergen administered. 32,33,63,64 Therefore, risk analysis for specific groups of patients (i.e., children vs. adults, comparison between various countries) had statistical limits which may influence results when comparing retrospective data from single studies or single protocol-specific datasets. 16,54 In conclusion, our study found similar milk EDs, but higher egg EDs compared with those previously found in other studies. Our data suggest that children with baked egg tolerance may be less exposed to the risk of severe reactions than children with baked milk tolerance.
They could be tolerating to egg traces and perhaps they could assume egg-PAL pre-packaged foods. Should our results be confirmed by studies comparing EDs differences between raw and baked food in pediatric population, the clinical attitude toward the two allergies could be differentiated. To confirm these data, we will evaluate if such difference in thresholds is present also in children allergic to baked egg/milk. For the moment, these data add to the knowledge about thresholds of reactivity for these foods, an essential characteristic of the hazard that allergens present to the food-allergic population. 65

Open access funding provided by Ospedale Pediatrico Bambino
Gesu.

AUTH O R CO NTR I B UTI O N S
RLV developed the concept, prepared the infrastructure of the study, elaborated the data, performed the statistical analysis, and contributed made the draft. AF initiated the concept and contributed made the first draft. CR, VC, LD, VF, MM, and AS administered the OFCs, evaluated the clinical data, included them in the database and contributed to the discussion. ALP collaborated to the inclusion of data in the database, collaborated to the statistical analysis, and contributed to the discussion. All authors reviewed and approved the final manuscript.

CO N FLI C T O F I NTE R E S T
All authors declare that they have no conflicts of interest to disclose. Note: ED01, ED05, ED10, ED50: cumulative amount of food protein predicted to cause an allergic reaction in 1%, 5%, 10%, or 50% of the food allergic population, respectively.