text:
First detection and evaluation of ivermectin resistance in Haemonchus contortus from sheep in Zhaosu, Xinjiang, China

Background: In 2017 in Zhaosu, Xinjiang, China, a series of sheep death suspected to have reduced the potency of ivermectin (IVM) have raised the vigilance of herdsmen. At the same time, veterinary practitioners warned of the potential risk of the lack of ecacy of IVM in the treatment of Haemonchus contortus. Methods: In this study, the ecacy of IVM against H. contortus was assessed under controlled experimental conditions by using larval developmental assay (LDA) and faecal egg count reduction test (FECRT) to determine whether the reduced ecacy in sheep in Zhaosu was caused by IVM resistance. Results: In total, 1560 adult worms were isolated from the abomasum of post mortem examination of dead sheep which received recent IVM (0.2 mg/kg) treatment before death and eggs were collected from female worms for usage in LDA. The LDA result revealed that signicant difference in resistance ratio (RR) was detected between the IVM-susceptible and -resistant isolates (6.09 for Zhaosu strain) of H. contortus. The Zhaosu strain was maintained in laboratory for FECRT, which showed that the ecacy of IVM was 0% at 1, 2 and 3 times dose, but at 4.8 times dose, it reached 87.5% at 14 days post-treatment. Conclusion: The ndings suggested that controlling H. contortus of small ruminants in Zhaosu requires an integrated parasite management strategy that should include monitoring of infection and anthelmintic resistance, and increasing the awareness of herdsmen on the importance of using the appropriate drug usage to improve nematode control in sheep. full ng/μl) larval assay larvae IC 50 values were 1.25 ng/μl and 0.205 ng/μl for Zhaosu and Haecon-5 strain, respectively. A comparison on the IC 50 value between Zhaosu strain Haecon-5 showed that resistance ratio (RR) is 6.09, Zhaosu used this resistance ratio as a reference and conducted FECRT using different Zhaosu strain.

Background
Zhaosu is located in the southwest of Yili Kazak Autonomous prefecture, Xinjiang, China, which is the only county without desert in Xinjiang. In Zhaosu, sheep are the most important small domestic ruminants which usually are reared in smallholder farming systems with great numbers. The local sheep industry has a variety of production systems including indoor and outdoor units with a sheep population of approximately 740000 (According to the statistics of Zhaosu animal husbandry and Veterinary Bureau in spring, 2020). After mid-May each year, most of the sheep breeding herd is kept outdoors and domesticated sheep often share pasture with other livestock species such as cattle and horses so that such environment signi cantly increases the chance of endo-parasite and ecto-parasite infections, and there is a potential risk for substantial clinical and economic losses caused by parasitism unless effective control is implemented.
Infection by Haemonchus contortus,a highly pathogenic and economically important parasite, is a serious problem for extensive small ruminant farming system in warm climatic regions on a worldwide basis [1]. Control of this parasite is achieved through a combination of anthelmintic treatment and management strategy, particularly rational administration of anthelmintic, good hygiene and land management to interrupt transmission [2,3]. Fortunately, three anthelmintic groups are available and have been used with great success for controlling this parasite in the past decades, including Benzimidazoles, Imidazothiazoles and Macrocyclic lactones [4,5,6]. Among them, Ivermectin (IVM) belonging to macrocyclic lactones is a highly effective, broad-spectrum and gold standard anthelminthic drug that has been used widely since 1990's to control both endo-and ecto-parasites in a wide variety of hosts [7]. Unfortunately, IVM resistance is becoming a serious problem to animal health and its resistance in H. contortus is a threat to sustainable livestock production as well as a growing worldwide concern [8,9]. Therefore, there is no guarantee that all chemicals will be uniformly effective in any region in the world.
Our scant understanding of the resistance mechanism of IVM in the H. contortus populations limits our ability to develop evidence-based mitigation strategies [10,11]. In spite of some reports on IVM resistance in H. contortus in other countries, there still lacks of information about IVM resistance in Xinjiang, China. However, concern of reduced IVM e cacy was raised to the Animal Husbandry and Veterinary Bureau in Zhaosu in recent years. Given the need for identifying the reasons causing several deaths of sheep in a herd in which IVM had been used regularly, the aim of our work was to con rm the presence of IVM resistance and to evaluate the resistance level in Zhaosu by using LDA in concert with FEC reduction tests (FECRT). The nal objective of this study was to provide a reliable data basis for local sheep industry.

Haemonchus contortus strains
Two isolates of H. contortus were used in the study:

Chemicals
The commercial IVM powder (Shyuanye Inc, Shanghai, China) was used as a source of IVM in LDA. The stock solution (1 mg/mL) was prepared in dimethyl sulfoxide (DMSO) and was diluted in DMSO to generate a series of working solutions. IVM Injection (Yuanzheng Pharmaceutical Co., Hebei, China) was used in the faecal egg count reduction test (FECRT).

Culture of larvae and morphological larval identi cation
The H. contortus female adults collected from the abomasum of dead sheep from Zhaosu were minced and then mixed with the sterilized sheep feces, which were then submitted to the College of Veterinary Medicine, Huazhong Agricultural University, Wuhan. The sample was incubated at 27 °C for 7 days, subsequently the third-stage larvae (L3) of Zhaosu strain of H. contortus were collected and morphologically identi ed. In order to further perform FECRT, around 6-months old goats free of parasites were infected with 8000 L3s per goat of Zhaosu strain and Haecon-5 strain, respectively. All experimental goats used in this study were housed in a cage without access to grass and feces.

Recovery of H. contortus eggs
Fresh H. contortus eggs were recovered from the faeces of each donor goat using a similar method from a previous report [12]. In short, approximately 100 g of faeces were macerated in a mortar with 800 ml of sugar water. The suspension was ltered using the 450 μm and 125 μm ne sieves and ltrate was covered with a piece of plastic lm by setting on the surface of suspension liquid for about 40 minutes.
Then, the plastic lm was removed and its surface was washed with normal saline. The collection solution was further ltered using the 50 μm and 38.5 μm ne sieves, the eggs on the surface of the sieves (38.5 μm) are collected and used immediately for larval development assays.

Larval development assay (LDA)
Laval development assay (LDA) was carried out as modi ed from the original method [13]. Brie y, the egg suspension (2 ml) recovered from faeces was added to T25 cell bottle containing 400 μl of growth medium (1 g yeast extract plus 90 ml dH 2 O; autoclaved for 15 min; added with 10 ml of 10× concentrated Earles's solution). Each bottle contains about 5000 eggs which were incubated at 27 °C for 24 h. On the next day, after the eggs hatched, 99 μl of larval suspension (~100 larvae) and 1 μl IVM working solution ( nal DMSO concentration in experimental group and control wells was 1% v/v with triplicate assay wells) were added to the wells of 96-well plates and incubated at 27 °C for another 6 d ( Fig received DMSO alone. At the end of the test, the number of L3s was counted in each well and the resistance ratio (RR) was calculated [14].
Faecal egg reduction test A faecal egg count reduction test (FECRT) was performed according to the World Association for the Advancement of Veterinary Parasitology (WAAVP) guidelines for evaluation of anthelmintic e cacy in ruminants [15]. In an attempt to further con rm the results of the LDA and the cause of death of sheep, four consecutive e cacy trials involving H. contortus experimental infections were carried out in goat. Thirty days after infection, goats in the treatment group were treated with IVM (IVM sub-cutaneous injection) at different doses including 0.2 mg/kg (standard dose), 0.4 mg/kg (2 times dose), 0.6 mg/kg (3 times dose) and 0.96 mg/kg (4.8 times dose), respectively. Goats in the control group were left untreated.
One day before the treatment and 14 d after the treatment, goats were checked for faecal egg counts (EPG) and the e cacy of the treatment was evaluated by the faecal egg count reduction test (FECRT), according to the formula, based on the guidelines recommended by WAAVP [15]. In addition, after 20 days of treatment at 0.4 mg/kg, the goats infected with Zhaosu strain were slaughtered to collect adults and calculate the worm burden. Percentage faecal egg count reduction (FECR) was determined according to the formula FECR (%) = 1-[T2/TI ] ×100 where T1 (before treatment) and T2 (post treatment) are the arithmetic mean egg counts in the group treated with IVM, respectively [16].

Case history
This investigation was carried out in July 12, 2017. According to the sheep owners, infected sheep became weak, lethargic and recumbent, despite treatment orally with IVM (0.2 mg/kg) before death. According to local veterinarians and farmers' statement, sheep in most farms in Zhaosu suffered from sporadic deaths in the different herd over the preceding two months (May to July 2017). Sheep were drenched orally with IVM (mainly) two to four times, each spring (April to May) and each autumn (September to November) in this region.

Post mortem examination
In July 2017, a 1-year-old male sheep which had been sick for several weeks and had received oral ivermectin (0.2 mg/kg) before its death was submitted to post mortem examination. Post mortem examination revealed that the faeces are rm and scant, the dead sheep was emaciated and malnourished, and had severe anemia. A total worm count got 1560 H. contortus worms in its gastric compartment (Fig. 2B). Sectional examination of the viscera didn't reveal other signi cant pathological changes and the cause of death was concluded to be severe anaemia and referable to the gastric H. contortus infestation.

Larval developmental assay
The results of dose response experiments with IVM against H. contortus larvae in vitro LDA revealed a good separation between susceptible Haecon-5 strain (control) and Zhaosu-resistant strain. In the control group, about 98% of the hatched larvae developed through to L3. In contrast, Zhaosu strain showed dose responses ranging from full development (0.2 ng/μl) to 100% inhibition on the development (12.5 ng/μl) with increasing IVM concentrations, indicating that the larval development assay was able to de ne the response of H. contortus larvae to the drug. The IC 50 values were 1.25 ng/μl and 0.205 ng/μl for Zhaosu strain and Haecon-5 strain, respectively. A comparison on the IC 50 value between Zhaosu strain and Haecon-5 strain showed that resistance ratio (RR) is 6.09, suggesting that Zhaosu strain was highly resistant to IVM. Then, we used this resistance ratio as a reference and conducted FECRT using different IVM doses to calculate IVM e cacy against Zhaosu strain.

Faecal egg count reduction test (FECRT)
In FECRT experiment, faecal eggs were detected in all ve groups including four experimental groups and one control group prior to administration of IVM, group mean egg counts still remain at similar level to that before treatment on different days after treatment at 600 ~ 1800 eggs per gram for the experimental and control groups (Table 1, Fig. 2C). Following IVM administration, in groups 1, 2 and 3, the FECs averaged EPG almost unchanged on day 14 post treatment (Table 1). However, for the 4th IVM administration group, average EPG ranged from 800 to 100 with the time going from Day 0 to Day 21 after treatment and FECR was 87.5% (Table 1). In this high dose group, although FEC signi cantly reduced, interestingly, it gradually increased with time and EPG reached 200 on day 28 post-treatment ( Table 1). Based on faecal egg count reduction, in four experimental groups, e cacies were 0%, 0%, 0% and 87.5% on day 14 post treatment, respectively (Table 1). In addition, on day 21 post-treatment in goat received 2 times dose IVM (0.4 mg/kg), 980 (46% vs 54% ) adult worms were isolated from the abomasum of post mortem examination (Fig. 2D).

Discussion
As H. contortus has shown a strong ability to develop resistance, IVM resistance in this parasite is extensively studied and well documented [17]. Today, most of the eld isolates of H. contortus from different continents were resistant to IVM, including North America [18], South America [19], Europe [20], Australasia [21], Africa [22] and Asia [23]. However, no registered cases of IVM resistance in H. contortus and IVM e cacy situation have been reported in Xinjiang. The present work constitutes the rst document of IVM resistance and follows reported treatment failures of IVM in sheep in Zhaosu, Xinjiang.
After mid April 2017, the practitioner noticed that there was a reduction in the average weight of sheep in different farms. They also observed that sheep suffered from exercise tolerance and subcutaneous oedema and sporadic death across different farms, and mortality of sheep over two months was about 10%. The clinical symptoms, morbidity and mortality of sheep over two months were suspected to be caused by haemonchosis, despite long-term regular oral drenching with IVM. The history of IVM usage and the post mortem ndings strongly suggested IVM resistance had led to a failure of H. contortus control in this region. Post mortem ndings and experimental results from LDA and FECRT in the present study con rmed that this suspicion was true, and indicated that IVM was completely ineffective. Therefore, the high level of resistance to IVM in H. contortus detected in this study is deeply concerning.
In the present study, it is important to highlight that surviving females showed no signi cant suppression in fecundity at 1, 2 and 3 times dose, but reduction in egg generation has been observed from the rst day to 20 days after treatment at 4.8 times dose. Similar ndings (reduced fecundity by IVM treatment) were also obtained in sheep (Cooperia curticei) [24], cattle (Cooperia curticei) [25] and pig (Oesophagostomum dentatum) [16] nematodes. A recent study in the free-living nematode Caenorhabditis elegans showed that calcium signaling plays a crucial role in fertilization [26]. Interestingly, another study in larial human parasitic nematode Brugia malayi reported that some genes in the Ca 2+ signaling pathway and a number of other genes involved in fertilization were down-regulated following exposure to IVM, which caused adverse effects on ovulation [27]. In view of the above report, the reduction in ovulation of the female worm of H. contortus may also be related with dysregulation of genes involved in the fertilization after IVM exposure.
Development of IVM resistance in this case might have been affected by a number of factors, such as excessive and frequent treatment, under-dosing with IVM, and its management strategies [28,29]. Firstly, high frequency of treatment is long recognized as a major causal factor for anthelmintic resistance [30], especially when using the same drug for a long time. IVM had been used for at least the last seven years without the use of other classes of anthelmintics in Zhaosu. There is no doubt that many livestock owners routinely use IVM which may remove su cient nematodes to achieve a good clinical result, but also exacerbate the resistance situation by allowing resistant worms to survive. We suggested that it could be the main reason for the prevalence of IVM-resistant worms in this region which would steadily increase over time. Secondly, lack of farmers' awareness to the importance of ensuring appropriate doses and an underestimation of animal weights that appears to have been common in Zhaosu, therefore, inadequate dose was a potential factor for resistance to develop [31,32,33]. Thirdly, the absence of detection and managing IVM resistance in the long term caused the selection pressure for the resistance in H. contortus to increase with time. In the majority of cases, FECRT is the only currently available method of simultaneously assessing the e cacy of anthelmintics and resistance level of species, especially in remote pastoral areas like Zhaosu. However, lack of technical support and epidemiological survey appears to be a barrier to its wide adoption and livestock owners are not aware of the anthelmintic resistance situation in their properties. Based on the results of this study, more comprehensive investigation are needed in the future to further con rm the clinical and economic impact of H. contortus IVM resistance on sheep in the local production system.

Conclusion
In conclusion, this study is the rst report on IVM resistance in H. contortus in Xinjiang, which is the initial step and important opportunity to understand and discuss the wider issue of anthelmintic usage and e cacy on sheep farms in Zhaosu. Local veterinarians and farmers should be aware of the relationship between the infection rate and the severity of haemonchosis and the possibility of resistance development with long term application of IVM in sheep production systems in order to implement more cost-effective control strategy. Therefore, wider investigation based on worm egg counts for worm prevalence, and faecal egg count reduction testing for drug resistance in the eld in sheep is timely in order to detect the emergency of drug resistance as early as possible and mitigate possible risk factors for further anthelmintic resistance development.  Figure 1 Map of Xinjiang, China: Red star represents exact location of the Zhaosu County which H. contortus specimens were collected. Note: The designations employed and the presentation of the material on this map do not imply the expression of any opinion whatsoever on the part of Research Square concerning the legal status of any country, territory, city or area or of its authorities, or concerning the delimitation of its frontiers or boundaries. This map has been provided by the authors.


extfieldsofstudy: Medicine Biology

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

text:
Genomic Prediction of Kernel Zinc Concentration in Multiple Maize Populations Using Genotyping-by-Sequencing and Repeat Amplification Sequencing Markers

Enriching of kernel zinc (Zn) concentration in maize is one of the most effective ways to solve the problem of Zn deficiency in low and middle income countries where maize is the major staple food, and 17% of the global population is affected with Zn deficiency. Genomic selection (GS) has shown to be an effective approach to accelerate genetic gains in plant breeding. In the present study, an association-mapping panel and two maize double-haploid (DH) populations, both genotyped with genotyping-by-sequencing (GBS) and repeat amplification sequencing (rAmpSeq) markers, were used to estimate the genomic prediction accuracy of kernel Zn concentration in maize. Results showed that the prediction accuracy of two DH populations was higher than that of the association mapping population using the same set of markers. The prediction accuracy estimated with the GBS markers was significantly higher than that estimated with the rAmpSeq markers in the same population. The maximum prediction accuracy with minimum standard error was observed when half of the genotypes were included in the training set and 3,000 and 500 markers were used for prediction in the association mapping panel and the DH populations, respectively. Appropriate levels of minor allele frequency and missing rate should be considered and selected to achieve good prediction accuracy and reduce the computation burden by balancing the number of markers and marker quality. Training set development with broad phenotypic variation is possible to improve prediction accuracy. The transferability of the GS models across populations was assessed, the prediction accuracies in a few pairwise populations were above or close to 0.20, which indicates the prediction accuracies across years and populations have to be assessed in a larger breeding dataset with closer relationship between the training and prediction sets in further studies. GS outperformed MAS (marker-assisted-selection) on predicting the kernel Zn concentration in maize, the decision of a breeding strategy to implement GS individually or to implement MAS and GS stepwise for improving kernel Zn concentration in maize requires further research. Results of this study provide valuable information for understanding how to implement GS for improving kernel Zn concentration in maize.

Enriching of kernel zinc (Zn) concentration in maize is one of the most effective ways to solve the problem of Zn deficiency in low and middle income countries where maize is the major staple food, and 17% of the global population is affected with Zn deficiency. Genomic selection (GS) has shown to be an effective approach to accelerate genetic gains in plant breeding. In the present study, an association-mapping panel and two maize double-haploid (DH) populations, both genotyped with genotypingby-sequencing (GBS) and repeat amplification sequencing (rAmpSeq) markers, were used to estimate the genomic prediction accuracy of kernel Zn concentration in maize. Results showed that the prediction accuracy of two DH populations was higher than that of the association mapping population using the same set of markers. The prediction accuracy estimated with the GBS markers was significantly higher than that estimated with the rAmpSeq markers in the same population. The maximum prediction accuracy with minimum standard error was observed when half of the genotypes were included in the training set and 3,000 and 500 markers were used for prediction in the association mapping panel and the DH populations, respectively. Appropriate levels of minor allele frequency and missing rate should be considered and selected to achieve good prediction accuracy and reduce the computation burden by balancing the number of markers and marker quality. Training set development with broad phenotypic variation is possible to improve prediction accuracy. The transferability of the GS models across populations was assessed, the prediction accuracies in a few pairwise populations were above or close to 0.20, which indicates the prediction accuracies across years and populations have to be assessed in a larger breeding dataset with closer relationship

INTRODUCTION
Known as "hidden hunger, " micronutrient malnutrition is mainly prevalent among pregnant women and infants in the low and middle income countries (LMIC), where people rely mostly on cereal-based diets (Diepenbrock and Gore, 2015). Low levels of micronutrients, including zinc (Zn), iron and provitamin A, lead to malnutrition-related health impairments (Cakmak, 2002;Tiwari et al., 2015). According to the World Health Organization, Zn deficiency affected 17% of the global population 1 . Zn micronutrient deficiency, prevalent among young children in developing countries, is associated with decreased immune-competence and increased rates of infectious diseases, which have been reported as an extensive foodrelated primary health problem in LMIC (Gibson, 1994;White and Broadley, 2009). Biofortification is a promising approach to improve micronutrient malnutrition through breeding and biotechnology, and enrich the micronutrient in the food by develop new varieties (Bouis and Saltzman, 2017).
The HarvestPlus project, a CGIAR research program, has been working to micronutrient malnutrition through bio-fortification of staple crops (Ortiz-Monasterio et al., 2007). Maize is one of the target crops of HarvestPlus, and the most important staple food for millions of people in major developing countries in sub-Saharan Africa, Latin America, and Asia (Jones and Thornton, 2003;Beyene et al., 2016). In maize, the baseline of kernel Zn concentration is about 20 mg/kg, the breeding target established by HarvestPlus project was 33 mg/kg, assuming the estimated average requirement of 1,860 µg per day of Zn in maize (Bouis and Welch, 2010). Therefore, an increase of at least 13 mg/kg is targeted by breeding. This target is achievable, due to the significant genetic variation for kernel Zn concentration exists in tropical maize germplasm, ranging from 4 to 96 mg/kg (Bänziger and Long, 2000;Ortiz-Monasterio et al., 2007;Prasanna et al., 2011Prasanna et al., , 2020Hindu et al., 2018). Enriching the kernel Zn concentration in maize through bio-fortification is one of the most effective ways to solve the problem of Zn deficiency for pregnant women and young children living in the above-mentioned areas. Dissecting the genetic architecture of kernel Zn concentration in maize with genome-wide molecular markers will allow breeders to improve their breeding efficiency by facilitating the introgression of the related genes into low Zn germplasm through marker-assisted selection or genomic selection (GS). Several studies have been conducted to dissect 1 www.harvestplus.org the genetic architecture of kernel Zn concentration in maize (Qin et al., 2012;Simic et al., 2012;Baxter et al., 2013;Jin et al., 2013;Hindu et al., 2018). Jin et al. (2013) identified five significant QTL and ten meta-QTL in 218 F 2:3 maize families. Hindu et al. (2018) detected 20 SNPs significantly associated with kernel Zn concentration in maize by implementing association mapping in a collection of 923 inbred lines, and 11 of those SNPs were validated in 3 DH populations by single marker linkage mapping analysis.
Genomic selection has been shown to be an effective approach to accelerate genetic gains in maize breeding (Lorenzana and Bernardo, 2009;Lian et al., 2014;Zhang et al., 2015Zhang et al., , 2017Cao et al., 2017). Highly variable prediction accuracy levels have been reported in plants depending on the training population size, the relationship between the training and the prediction sets, trait complexities, marker densities, and genotyping platforms (Zhao et al., 2012;Bian and Holland, 2017;Norman et al., 2018). Several studies have shown similar results, i.e., that prediction accuracy increases as trait heritability, size of training set, and number of markers in various types of maize populations increase (Heffner et al., 2009;Hickey et al., 2014;Liu et al., 2018;Norman et al., 2018). Norman et al. (2018) showed that prediction accuracy can be improved by broadening the genetic diversity within the training set, particularly when relatedness between training and validation sets is low. Genomic prediction analyses were conducted on a collection of 284 maize inbred lines, which were genotyped with both 1,148 and 55,000 SNPs. Results indicated that the prediction accuracies increased as the number of markers used across all the trait-environment combinations increased (Crossa et al., 2010;Gonzalez-Camacho et al., 2012).
An economical genotyping platform is always required in order to make GS more cost-effective. GBS, a next-generation sequencing technology, is a high-throughput, multiplex and short-read sequencing approach that reduces genome complexity via restriction enzymes and generates high-density genomewide markers (∼1 million) at a low cost per sample by tagging randomly shared DNA fragments from different samples with unique, short DNA sequences (barcodes) and pooling samples into a single sequencing channel (Elshire et al., 2011;Wu et al., 2016). Several studies have indicated that GBS is a promising genotyping platform for GS applications (Poland et al., 2012;Crossa et al., 2013;Zhang et al., 2015;Yu et al., 2016). Zhang et al. (2015) showed that GBS outperformed lowdensity SNPs for both complex and simple traits evaluated under stress conditions with low-to-moderate heritability in 19 tropical maize bi-parental populations evaluated in multi-environment trials. Rio et al. (2019) obtained good prediction accuracies for grain moisture, grain yield, yield index and male flowering in a collection of 389 dent kernel type maize inbred lines that were genotyped with GBS and other sequencing technology. Developing a high-throughput genotyping platform with high quality, flexibility, and affordable genotyping cost is still critical for implementing GS routinely in large-scale breeding programs. rAmpSeq is a simple, robust platform for designing primers, PCR amplification, and high-throughput multiplex sequencing, which allows hundreds to thousands of markers to be scored for less than $5 per sample (Buckler et al., 2016). rAmpSeq is specifically tailored to GS approaches (Voss-Fels et al., 2019). However, GS applications that include genotyping the training and validation sets with rAmpSeq markers have not been reported until now.
In the present study, an association-mapping panel and two maize DH populations, genotyped with GBS and rAmpSeq markers, were used to estimate the genomic prediction accuracy of kernel Zn concentration in maize. The main objectives were to: (1) estimate the genomic prediction accuracy for kernel Zn concentration in different maize populations using GBS and rAmpSeq markers; (2) compare the genomic prediction accuracies of kernel Zn concentration within and across multiple maize populations estimated by the different genotyping platforms; (3) assess the effect of training population size (TPS), marker density (MD) and marker quality on genomic prediction accuracy estimation; and (4) explore training population development base on the phenotypic variation of the target trait.

Plant Materials
An association-mapping panel and two DH populations were used to perform genomic prediction analyses in the current study. The association-mapping panel, designated Drought Tolerant Maize for Africa (DTMA), consists of 300 lowland tropical and mid-altitude tropical inbred lines. These lines originated from different CIMMYT maize breeding programs and have abundant genetic variation.
The two DH populations, designated DH1 and DH2, were derived from F1 crosses between two inbred lines. DH1 consisted of 108 lines, and the F1 cross was made between two elite maize inbred lines of CML503 and CLWN201. DH2 consisted of 143 lines, and the F1 cross was made between two elite maize inbred lines of CML465 and CML451. Plot size for all experiments was single row 2.5 m length, with 75 cm between rows, and 0.23 m between plants in each row. In each plot, six plants were self-pollinated, hand-harvested and hand-shelled to avoid any metal contamination, and bulked kernel samples from each plot were dried and sent for analysis at the maize nutritional quality analysis laboratory in Mexico. More details on the analyses of kernel Zn concentration have been previously described by Hindu et al. (2018).

Phenotypic Data Analysis and Heritability Estimation
In each of the three populations, least-squares means of genotypes were calculated across environments through the "lsmeans" function of the R program, version 5.6.1 (R Core Team, 2019). The least-squares means were proposed to analyze the data with unequal subclass numbers, and which are predictions from linear or mixed models (Lenth, 2016). For the locations with replicated data, mean values across the two replications were calculated for predictions of least-squares means. A linear mixed model was fitted to the data as follow: where Y ik is the mean performance of a certain genotype, µ is the overall mean effect, Env i is the effect of i th location, Gen k is the main effect of the k th genotype, ε ik is the error associated  with the i th location and the k th genotype, which is assumed to be normally and independently distributed, with mean zero and homoscedastic variance. All factors except genotype were set as random. Narrow-sense heritability was estimated as the ratio of additive genetic variance to total phenotypic variance: where V A was an estimate of the additive genetic variance, and V P was the total phenotypic variance. The total phenotypic variance was the sum of V A and V e , and V e was an estimate of the residual variance. The variance components of V A and V e were estimated based on the genomic relationship matrix in the rrBLUP package (Endelman, 2011) of the R program, version 5.6.1 (R Core Team, 2019). In both the DTMA panel and the two DH populations, the Pearson correlation coefficients of kernel Zn concentration between locations were estimated in each population using the R program version 5.6.1 as well.

Genotyping and Genotypic Data Analysis
For all the lines in each population used in the present study, leaf samples of each line were collected 3-4 weeks after seeding for DNA extraction with a CTAB procedure (CIMMYT Laboratory protocol, 2005). All the lines in the two DH populations and a FIGURE 2 | Distribution of MAF and missing rate across all the three populations before and after filtering: (A) MAF distribution of the GBS marker dataset before filtering; (B) MAF distribution of the GBS marker dataset after filtering; (C) MAF distribution of the rAmpSeq marker dataset before filtering; (D) MAF distribution of the rAmpSeq marker dataset after filtering; (E) missing rate distribution of the GBS marker dataset before filtering; (F) missing rate distribution of the GBS marker dataset after filtering.
subset of the DTMA panel of 236 inbred lines were sent to the Biotechnology Resource Center of Cornell University for both GBS and rAmpSeq.
A GBS protocol commonly used by the maize research community was applied in this study (Elshire et al., 2011). Genomic DNA was digested with the restriction enzyme ApeKI, and a DNA library was constructed in 96-plex and sequenced on Illumina HiSeq2000. Details of analyses of SNP calling and imputation have been previously described . Initially, 955,690 SNPs evenly distributed on maize chromosomes were called for each line; 955,120 of them were assigned to chromosomes 1-10, and 570 of them could not be anchored to any of the 10 maize chromosomes.
rAmpSeq is a simple, robust, and cost-effective genotyping strategy developed by Cornell University for large-scale GS projects. Details of the rAmpSeq primer pair information have been described by Buckler et al. (2016). DNA libraries were constructed in 3072-plex and sequenced on Illumina HiSeq2000, and each sequence tag was treated as a unique dominant marker. Initially, 7,595 dominant markers identified from the intergenic regions were called for in each of the genotyped samples.

Genomic Prediction Analysis Within Each Population
Genomic prediction was performed in the rrBLUP package in the R program version 5.6.1 (Endelman, 2011). The mixed model is described as: where y is the vector (n × 1) of observations, X is the vector (n × 1) of individuals and β is the fixed effects, ε is the vector (n × 1) of independently random residuals with assumed distribution N (0, Iσ ε 2 ), Z is the design matrix (n × m) for random effects, and u is the vector of random effects with u ∼ N (0, Kσ u 2 ), K being an identity matrix in this case (Endelman, 2011;Liu et al., 2018). In addition, n is the number of individuals, and m is the number of markers.
A fivefold cross-validation scheme was used to generate the training and validation sets and assess the prediction accuracy within each population, where the average value of the correlations between the true breeding value and the genomic estimated breeding values was defined as genomic prediction accuracy. In each population, the data were divided into two subsets, with 80% of the lines randomly selected and assigned to the training set and the remaining 20% assigned to the validation set. In total, 150 replications of cross-validation were performed for each population. In the GBS dataset, TASSEL version 5.0 (Bradbury et al., 2007) was used to filter the markers in each of the three populations with a minor allele frequency (MAF) greater than 0.05 and a missing rate less than 20%. In total, 262,919, 65,430, and 46,426 GBS markers were selected for further genomic prediction analyses in the DTMA panel, DH1 population, and DH2 population, respectively. In the rAmpSeq dataset, markers with an MAF greater than 0.05 were filtered in each of the three populations. In total, 6,150, 3,859, and 2,795 markers were selected for further genomic prediction analyses in the DTMA panel, DH1 population, and DH2 population, respectively. In the same population, the prediction accuracies estimated from the GBS and rAmpSeq datasets were compared, and a t-test was conducted to obtain the significance.

Effect of Training Population Size (TPS), Marker Density (MD) and Marker Quality on Genomic Prediction Accuracy Estimation
To assess the effect of TPS on the estimation of genomic prediction accuracy, the training population was set from 10 to 90%, with an interval of 10%, in each of the three populations. The number of markers used for prediction in each population were same as those used in the cross-validation analyses. In total, the analysis was repeated 100 times in each population.
To evaluate the effect of MD on the estimation of genomic prediction accuracy, the number of markers varying from 10 to all markers (i.e., 10, 50, 100, 300, 500, 1,000, 3,000, 5,000, 10,000, 50,000, and all markers) were used for genomic prediction analyses in each of the three populations. In the GBS datasets, all markers with an MAF greater than 0.05 and a missing rate less than 20% were filtered. In the rAmpSeq datasets, all markers with an MAF greater than 0.05 were filtered. The fivefold crossvalidation scheme was repeated 100 times in each population with different marker datasets.
To examine the effect of marker quality on the estimation of genomic prediction accuracy, different levels of MAF and missing rates were used to filter the marker datasets and control marker quality. In GBS datasets, markers were filtered with the combinations between MAF and missing rate in each population; MAF ranging from 0.10 to 0.40, with an interval of 0.10; missing rates ranging from 0 to 80%, with an interval of 20%. In rAmpSeq datasets, markers were filtered with MAF in each population, and MAF ranging from 0.10 to 0.40, with an interval of 0.10. The fivefold cross-validation scheme was used to compare the prediction accuracies estimated from the marker datasets with different quality levels. In total, the analysis was repeated 150 times in each population.

Training Set Development Based on the Phenotypic Variation of the Target Trait
Training sets were formed according to the phenotypic variation information of the target trait. Five scenarios were simulated and compared in each of the three populations, where the training set was formed by sampling the same percentage of genotypes with random selection (Random), with selection from the top tail (Top), with selection from the bottom tail (Bottom), with selection from the middle part (Middle), and with selection from the two tails (Two tails). In each scenario, the validation set was the whole population, and the training set ranged from 10 to 90%, with an interval of 20%. In each of the three populations, a total of 25 combinations and comparisons were conducted between the five scenarios and the five percentage levels of the training set.

Genomic Prediction Analysis Between Pairwise Populations
Among the three populations, genomic prediction analyses between pairwise populations were performed, when one population was used as training set to predict the other population as validation set, the correlations between the true breeding value and the genomic estimated breeding values of the validation set was defined as genomic prediction accuracy. In the GBS dataset, TASSEL version 5.0 (Bradbury et al., 2007) was used to filter the markers across the three populations with a MAF greater than 0.05 and a missing rate less than 20%. In the rAmpSeq dataset, markers with an MAF greater than 0.05 were filtered across the three populations. The genomic prediction analyses between pairwise populations were performed with 137,593 and 6,005 markers selected from the datasets of GBS and rAmpSeq, respectively.

Genomic Prediction Analysis With the Significantly Associated Markers
Genomic prediction analyses with significantly associated markers were performed to simulate marker assisted selection (MAS). In the previous study of Hindu et al. (2018), 20 GBS SNPs significantly associated with kernel Zn concentration in maize were identified in a collection of 923 inbred lines, and the DTMA panel used in the present study is a subset of these 923 inbred lines. In total, 11 of these 20 significantly associated SNPs were validated in three DH populations, and the two DH populations used in the present study are a subset of these three DH populations. A fivefold cross-validation scheme was used to assess the prediction accuracy of MAS within each population, when the 11 validated GBS SNPs were selected to perform genomic prediction. The comparison between GS and MAS was only applied in the GBS dataset, which is not applied in the rAmpSeq dataset, because of the lack of information about the physical position of the rAmpSeq markers.

Phenotypic Variation
Kernel Zn concentration in the all genotypes across the three populations ranged from 16.87 to 39.53 mg/kg, with an average value of 26.11 mg/kg ( Table 1). The average value of kernel Zn concentration was 27.11, 24.59, and 25.59 mg/kg in the DTMA panel, DH1 population, and DH2 population, respectively. The DTMA panel had the widest range of variation among the three populations, although the standard deviation within the two DH populations was similar to or higher than that within the DTMA panel. In the two DH populations, the greater standard deviation values indicated the larger genotype-by-environment interactions in a mixture of replicated and un-replicated trials, and the lower correlations between locations. In total, 2.87% of 487 inbred lines (i.e., 14 lines) reached 33.00 mg/kg (Figure 1), the target level of HarvestPlus, through bio-fortification (Bouis and Welch, 2010). Among the fourteen lines, nine were from the DTMA panel, four were from DH2, and only one was from DH1. The estimated narrow-sense heritabilities of the three populations were moderate to high, and the highest heritability (0.84) was observed in the DTMA panel, while the lowest heritability (0.62) was observed in DH2. In each population, the Pearson correlation coefficients of kernel Zn concentration between locations were moderate, ranging from 0.43 to 0.62 (details not shown).

Distribution of MAF and Missing Rate
In the GBS dataset across the three populations, the average MAF was 0.09, with continuous distribution classes from 0 to 0.50 at intervals of 0.05. Before filtering, 64.80% of the markers had an MAF < 0.05. In the other intervals (0.05-0.50), the percentages of markers in each interval were below 10%, ranging from 2.70 to 9.20%. After filtering, average MAF across three population had a significant increase to 0.23. In the GBS dataset across the three populations, the average missing rate across all the markers was 0.29 before filtering, and it decreased to 0.08 after filtering (Figure 2).
In the rAmpSeq dataset across the three populations, the average MAF was 0.19 and 35.40% of the total markers had an MAF < 0.05. In the other intervals (0.05-0.50), the percentages of markers in each interval ranged from 11.10 to 22.80% (Figure 2). After filtering, MAF across three population increased to 0.23.

Genomic Prediction Accuracies Estimated From the Fivefold Cross-Validation Schemes Within Each Population
Genomic prediction accuracies estimated from the fivefold crossvalidation schemes for all three populations are shown in Figure 3, where the prediction accuracies were moderate and varied across populations and genotyping platforms. Among the three populations, the lowest prediction accuracy was observed in the DTMA panel across both the GBS and rAmpSeq marker datasets. In the same population, the prediction accuracy estimated from the GBS marker dataset was higher than that estimated from the rAmpSeq marker dataset, and the difference was significantly. The prediction accuracies estimated from the GBS marker dataset were 0.40, 0.64, and 0.65 for the DTMA panel, DH1 population, and DH2 population, respectively. The prediction accuracies estimated from the rAmpSeq marker dataset were 0.35, 0.61, and 0.50 for the DTMA panel, DH1 population, and DH2 population, respectively (Figure 3). In the same population, the prediction accuracies estimated from the GBS marker dataset had less variation than those estimated from the rAmpSeq marker dataset.

Effect of TPS, MD, and Marker Quality on the Estimation of Genomic Prediction Accuracy
In both the GBS and rAmpSeq datasets, the prediction accuracy increased continuously as the TPS increased across all the populations (Figure 4). In the GBS datasets, the prediction accuracy increased slightly in both the DTMA panel and the two DH populations, when the TPS increased from 50 to 90%. The smallest standard error was observed in prediction accuracy when 60% of the total genotypes were assigned as the training set in the DTMA panel. In the two DH populations, the smallest standard error was observed when 50% of the total genotypes were assigned as the training set. In the rAmpSeq datasets, a similar trend was observed in all the three populations, which indicated that 50-60% of the total genotypes assigned as the training set can achieve good prediction accuracy. The effect of marker density on the estimation of prediction accuracy is presented in Figure 5. In the DTMA panel, the prediction accuracy continuously increased as the number of markers increased across both the GBS and rAmpSeq datasets. In the DH populations, a slight increase was observed in prediction accuracy when the number of markers increased from 300, the prediction accuracies nearly reached a plateau at 500 markers in both the GBS and rAmpSeq datasets. The results indicated that a larger number of markers is required to obtain higher prediction accuracy in populations with greater genetic diversity.
The result of prediction accuracies estimated in all the populations under the different marker datasets filtered with combinations of MAF and missing rate, or only with MAF, is presented in Tables 2, 3. This result showed a slight difference when compared with the prediction accuracies estimated with an MAF greater than 0.05 and a missing rate less than 20% in fivefold cross-validation schemes. In the GBS datasets across all the populations, MAF had a greater effect than the missing rate on the estimation of prediction accuracy, especially an MAF interval of 0.40-0.50 in the DH populations. When the MAF interval was 0.40-0.50, the average prediction accuracy estimated across all levels of missing rate was 0.40, 0.49, and 0.56 in the DTMA panel, DH1 population, and DH2 population, respectively, while the average prediction accuracy estimated across all missing rate levels in other MAF intervals ranging from 0.00 to 0.40, was 0.41, 0.65, and 0.66 in the DTMA panel, DH1 population, and DH2 population, respectively. Similar trends were observed in the rAmpSeq dataset across all the populations, a decreases in the prediction accuracy was shown in the MAF interval of 0.40-0.50. Missing rate had a minor effect on the estimation of prediction accuracy, and a very slight difference in prediction accuracy was observed among the different levels of missing rate, except for the DTMA panel with a missing rate of 0%. Across all the MAF intervals, the average prediction accuracy estimated at the missing rate level of 0% was 0.36 in the DTMA panel, while the average prediction accuracy estimated in the DTMA across the other levels of missing rate was 0.42, with a range of 0.39 to 0.43.

Training Set Development Based on the Phenotypic Variation of the Target Trait
For all three populations, the result of prediction accuracies estimated in the 25 combinations between the five scenarios and the five training set percentages are presented in Figure 6 for both the GBS and rAmpSeq datasets. Across all five scenarios and marker datasets, the prediction accuracy increased in all the populations as the increase of the percentages of training set. For example, the prediction accuracy estimated with the GBS dataset in the "Top" scenario was 0.04, 0.29, 0.58, 0.69, and 0.86, when the training set percentage in the DTMA panel was 10, 30, 50, 70, and 90%. Across all the training set percentages in all the populations, the "Two tails" scenario outperformed the other four scenarios in both the GBS and rAmpSeq datasets. For example, the prediction accuracy in the DTMA panel estimated with the GBS dataset at a training set percentage of 50% was 0.91, 0.75, 0.55, 0.30, and 0.44 for the Two tails, Random, Top, Middle, and Bottom scenarios, respectively. Similar trends were also observed in the two DH populations. These results indicated that developing a training set with broad phenotypic variation is possible to improve the prediction accuracy.

Genomic Prediction Accuracies Estimated From the Pairwise Populations
Genomic prediction accuracies between pairwise populations are shown in Table 4, where the prediction accuracies were very low across all the pairwise populations and genotyping platforms. In the GBS dataset, the average prediction accuracy across the six pairwise populations was 0.04. The highest prediction accuracy value was 0.32, when the DTMA panel was used as training set to predict the DH1 population as validation set. The prediction accuracies in other five pairwise populations were close to zero. In the rAmpSeq dataset, the average prediction accuracy across the six pairwise populations was 0.08. The highest prediction accuracy value was 0.24, when the DH2 population was used as training set to predict the DH1 population as validation set. The second highest prediction accuracy value was 0.19, when the DTMA panel was used as training set to predict the DH1 population as validation set.

Prediction Accuracy of MAS Estimated With the Significantly Associated SNPs
Genomic prediction accuracies estimated with the significantly associated SNPs are shown in Table 5

DISCUSSION
The main advantage of GS over phenotype-based selection is that it can accelerate the genetic gain per unit time and unit cost by reducing the selection cycle time and the phenotyping cost. However, the prediction accuracy must be high enough for GS to be effective. In the present study, an association-mapping panel and two maize DH populations genotyped with GBS and rAmpSeq markers were used to estimate the genomic prediction accuracy of kernel Zn concentration in maize. Results indicated that the prediction accuracies of kernel Zn concentration in maize were moderate to high and varied across populations and genotyping platforms. The prediction accuracy of kernel Zn concentration in the association panel estimated with GBS and rAmpSeq markers was 0.40 and 0.34, respectively. In the two DH populations, the prediction accuracies of kernel Zn concentration estimated with GBS markers ranged from 0.64 to 0.65, while the prediction accuracies of kernel Zn concentration estimated with rAmpSeq markers ranged from 0.50 to 0.61.
In the same population, the prediction accuracy estimated from the GBS marker dataset was higher than that estimated from the rAmpSeq marker dataset, and the difference was significantly. Several genomic prediction studies using GBS markers were implemented successfully in maize to improve various traits with different levels of genetic complexity, where the genotyping cost was at least $35 per sample at the 96plex level, or $13 per sample at the 384-plex level . However, a high-throughput genotyping platform with affordable genotyping cost is still required for implementing GS routinely in the large-scale breeding programs. A low-cost genotyping platform allows GS to be more cost-effective, and makes it feasible to replace expensive phenotyping with cheaper genotyping. The total breeding population size increases under the same budget by phenotyping a lower number of breeding lines in the training set and genotyping a greater number of breeding lines in the prediction set. Therefore, the selection intensity increases to accelerate the genetic gain per unit cost. rAmpSeq is specifically tailored to GS approaches, with a cost of $5 per sample. To the best of our knowledge, this is the first report of a genomic prediction study in maize using rAmpSeq markers, and the results of this study showed that GS using rAmpSeq markers is a cost-effective approach to improve the kernel Zn concentration in maize through bio-fortification. GS using rAmpSeq markers is also being implemented in CIMMYT maize breeding programs for improving grain yield and kernel Zn concentration simultaneously. The preliminary cost-benefit analysis showed that a breeding strategy that implements phenotype-based selection and GS stepwise could reduce the breeding cost up to 50% compared with phenotype-based selection only, equivalent to double the total breeding population size under the same budget.
Results of this study showed that the prediction accuracies continuously increased as the TPS increased in all the populations. Across the two genotyping platforms, relatively high prediction accuracies with the smallest standard error were observed in all populations, when 50 to 60% of the total genotypes were used as a training set. Results of this study are consistent with previous reports Cao et al., 2017), where the results also recommend phenotyping and genotyping as few as half of the genotypes as the training set to achieve good prediction accuracy of the target trait. In addition, prediction accuracy could be further improved by polling multiple populations as the training set to increase the population size, which will be assessed in further studies.
The transferability of the GS models across populations was assessed by estimating the prediction accuracies between pairwise populations, when one population was used as training set to predict the other population as validation set. The prediction accuracies were very low across all the pairwise populations and genotyping platforms, the prediction accuracies in a few pairwise populations are above or close to 0.20, which indicates that the genomic prediction accuracies estimated across populations could be improved furtherly by increasing the TPS and strengthening the relationship between the training and prediction sets. The transferability of the GS models across populations and across years are being tested with a larger dataset from a maize breeding program, where the TPS is bigger and the relationship between the training and prediction sets is closer. The preliminary results are promising, and more details will be reported in further studies.
How marker density affects prediction accuracy has been investigated in several previous studies, and the number of markers required to achieve good prediction accuracy may vary depending on the extent of linkage disequilibrium between markers and QTL, population types, and the genetic complexity of the target trait. In this study, the results showed that the prediction accuracy of kernel Zn concentration in maize continuously increased as the number of markers increased across populations and genotyping platforms. Across the two genotyping platforms, the prediction accuracy reached a plateau at 3,000 markers in the DTMA panel, and at 500 markers in the two DH populations, which indicated that a larger number of markers are required to obtain higher genomic prediction in populations with greater genetic diversity. In the tropical maize association-mapping panel, the average linkage disequilibrium decay distance over all 10 chromosomes was less than 5 kb at r 2 = 0.1, and roughly an estimated 500,000 markers are required to ensure that at least one marker can be in linkage disequilibrium with trait-associated loci. The DH population has a clear genetic structure and finite chromosome recombination events (Smith et al., 2008). Therefore, less than 500 markers are enough to ensure that at least one marker can be in linkage with each gene-related locus.
There is a tradeoff between number of markers and marker quality, because marker quality becomes lower as the number of markers increase in a specific marker dataset. The number of markers affecting prediction accuracy has been investigated in several previous studies, but very few reports have been conducted to estimate the effect of marker quality on prediction accuracy estimation. Results of this study show that the prediction accuracy reached a plateau at 3,000 markers in the DTMA panel, and at 500 markers in the two DH populations. It indicated that more markers with lower quality have little effect on prediction accuracy improvement, and markers with lower quality could be the noise on improving prediction accuracy. Moreover, the computational burden increases when a higher number of markers are used for prediction. Appropriate levels of MAF and missing rate should be considered and selected to improve the prediction accuracy and reduce the computational burden by balancing the number of markers and marker quality.
Several previous studies estimated the genetic diversity of the training set with molecular markers, and assessed the genetic diversity of the training set on prediction accuracy improvement (Dos Santos et al., 2016;Norman et al., 2018). Results of this study indicated that the genetic diversity of the training set could be estimated not only with molecular markers, but also with phenotypic data. Therefore, prediction accuracy can be improved by developing a training set with broad phenotypic variation.
In the present study, GS outperformed MAS and showed higher genomic prediction accuracies within each population on predicting the kernel Zn concentration. These results showed that kernel Zn concentration in maize could be improved by implementing GS individually or by implementing MAS and GS in a stepwise fashion. As we discussed earlier, a breeding strategy that implements phenotype-based selection and GS stepwise could reduce the breeding cost up to 50% compared with phenotype-based selection only, equivalent to double the total breeding population size under the same budget, if the breeding target is to improve grain yield and kernel Zn concentration simultaneously. The long-term genetic gain of how much GS could provide to improve grain yield and increase the Zn concentration will be further assessed in CIMMYT maize breeding programs. Alternative, kernel Zn concentration in maize could be improved by implementing MAS and GS in a stepwise fashion, where MAS is implemented as forward breeding at an early generation on larger numbers of selection candidates, followed by GS at advanced stages of breeding on smaller number of selection candidates for further improvement. The decision of a breeding strategy to implement MAS and GS stepwise for improving kernel Zn concentration in maize requires further research.

DATA AVAILABILITY STATEMENT
The datasets generated in this study can be found in the following repository: https://data.cimmyt.org/dataset.xhtml?persistentId= hdl:11529/10548362.

AUTHOR CONTRIBUTIONS
XZ, LZ, and TD initiated and designed the overall study. NP-R, FS, and TD performed and coordinated the field experiments and phenotypic data collection. XZ, MO, ML, and BP contributed to the genotypic data generation. RG, EM, DY, YR, AZ, and JC carried out the statistical analysis. RG, XZ, LZ, and TD interpreted the results and wrote the manuscript. All authors contributed to manuscript editing.


extfieldsofstudy: Medicine Chemistry

text:
Assessment of Volatiles and Polyphenol Content, Physicochemical Parameters and Antioxidant Activity in Beers with Dotted Hawthorn (Crataegus punctata)

Beer with the addition of dotted hawthorn (Crataegus punctata) fruit and juice was prepared and analysed. The content of carbohydrates, glycerol and ethanol in beers was determined by high-performance liquid chromatography (HPLC). Analysis of the total content of polyphenols was also performed using the Folin-Ciocalteu method, as well as determining antioxidant capacity by DPPH• and ABTS+• assay, and the ability to reduce iron ions by FRAP assay. Content of volatile compounds was analysed by means of solid-phase microextraction and gas chromatography coupled with mass spectroscopy. Beers with addition of hawthorn, both juice and fruit, had higher antioxidative potential and higher polyphenols concentration compared to control beer. The content of polyphenols in beers was in the range 200.5–410.0 mg GAE/L, and the antioxidant activity was in the range of 0.936–2.04 mmol TE/L (ABTS+• assay), 0.352–2.175 mmol TE/L (DPPH• assay) and 0.512–1.35 mmol TE/L (FRAP assay). A sensory evaluation of beers was also carried out. Beer with hawthorn fruit addition obtained the best scores in sensory analysis for criteria such as aroma, taste and overall quality. This beer had the highest content of volatile compounds (287.9 µg/100 mL of beer), while the control beer had lowest concentrations (35.9 µg/100 mL of beer).

Introduction
Beer polyphenols come from malt (70-80%) and hops (20-30%) [1]. Polyphenols and melanoidins affect the antioxidant activity of beer, as well as its sensory properties [2]. Antioxidant compounds play a role in reducing the amount of free radicals present in human body. The imbalance between antioxidants and pro-oxidants is called oxidative stress, which is a cause of many diseases. Exogenous antioxidants, introduced into the body with food, can scavenge free radicals, chelate metal ions, or inhibit pro-oxidative enzymes, protecting our body from the negative effects of oxidative stress [3].
Moderate consumption of beer has a positive effect on some biomarkers associated with human health [4,5]. Food rich in polyphenols can reduce the possibility of developing many civilization-related diseases such as arteriosclerosis, cancer and neurodegenerative diseases. This means that there is a public interest in products with high content of polyphenols [6,7]. Quifer-Rada et al. identified

Research Material
The research material was American Saison-style beer in three variants: without any additions (BC), enriched with the addition of hawthorn fruit (BF) and enriched with hawthorn fruit juice (BJ). The raw materials were Pilsen malt, wheat malt, Vienna malt, Caramel Pale malt, Light Caramel malt, acidified malt. All malts were acquired from the Viking Malt company (Strzegom, Poland). Amarillo and Calypso hops (USA) were purchased from Twój Browar company (Wrocław, Poland). The hawthorn fruit (Crataegus punctata) of the variety "Aurea" and juice pressed from this fruit were added in the production process. The juice was obtained from fruits by the use of a manual winepress.

Brewing Technology
Mashing was carried out under laboratory conditions with an infusion system in the following conditions: 64 • C for 60 min and 72 • C for 15 min. Next, the whole mash was heated to 78 • C to inactivate the malt enzymes, filtered and 25 L wort was obtained. The wort was boiled for 75 min with addition of hops in three doses: first, Calypso (15 g), 75 min of boiling; second, Amarillo (15 g), 30 min of boiling; third, Amarillo (15 g) and Calypso (10 g), 15 min of boiling. The hopped wort was cooled to 25 • C, filtered and aerated. The initial extract content was set at 12.5% (w/w) measured using Densito 30PX densimeter (Mettler Toledo, DC, USA). Fermentation was carried out in fermentation vessels at 25 • C. Hawthorn juice or fruits were added to beers at a dose of 10% (w/v) in the last 3 days of fermentation. The resulting beer was aged in bottles for 3 weeks.

Physicochemical Parameters
The pH value of beer was measured using Mettler-Toledo MP 220 (Greifensee, Switzerland). The extract content was measured at 20 • C with a Densito 30P density meter (Mettler-Toledo, Columbus, OH, USA). The measurements were performed in triplicate.

Determination of Carbohydrate Profile, Ethanol and Glycerol Content
The sugar profile and the content of ethanol and glycerol were examined by means of High-Performance Liquid Chromatography (HPLC) [26]. Degassed and centrifuged (5500 rpm, 10 min) beer samples were diluted two-fold (1:1) with ultrapure water and filtered through nylon syringe filters with 0.45 µm pore size for chromatographic vials. The samples were then analysed using a Prominence liquid chromatography system (Shimadzu Corp., Kyoto, Japan) equipped with a Rezed ROA-Organic Acid H + column (300 × 4.6 mm) from Phenomex (Torrance, CA, USA). The following measurement parameters were used: sample volume: 20 µL, separation temperature: 60 • C, mobile phase flow rate: 0.6 mL/min, mobile phase: 0.005 M H 2 SO 4 , detection temperature: 50 • C. The concentration of ethanol, glycerol, dextrins, maltose, glucose and maltotriose was based on five-point calibration curves using Chromax 10.0 software (Pol-Lab, Wilkowice, Poland). All measurements were performed in triplicate.

Analysis of Total Polyphenol Content
The total polyphenol content of the beer was determined using the Folin-Ciocalteu (F-C) spectrophotometric method [27]. Degassed and centrifuged beer samples or fruit juice samples (0.1 mL) and 0.2 mL of F-C reagent were pipetted into cuvettes. After 3 min, 1 mL of a 20% aqueous solution of sodium carbonate (Na 2 CO 3 ) and 2 mL of distilled water were added. The absorbance at 765 nm was measured after 1 h using Beckmann DU650 spectrophotometer (Brea, CA, USA) and the results were expressed as mg of gallic acid equivalents (GAE) per L of beer. Data were expressed as the mean value for three measurements. The antiradical activity was determined using a DPPH • radical assay [28]. Samples of beer or fruit juice (0.1 mL) were mixed with 2 mL of 0.04 mmol/L DPPH • in ethanol and 0.4 mL of distilled water. After 10 min of incubation at room temperature, the absorbance was measured with a spectrophotometer at 517 nm using disposable polystyrene cuvettes. A calibration curve was prepared with Trolox solution (0.005 mmol/L). The data were expressed as Trolox equivalent (TE) of antioxidative capacity per 1 L of the beer (mmol TE/L). All measurements were performed in triplicate.

Ferric Reducing/Antioxidant Power (FRAP) Assay
The FRAP assay is based on the reduction of ferric 2,4,6-tris(2-pyridyl)-1,3,5-triazine [Fe (III)-TPTZ] to the ferrous complex at low pH, followed by spectrophotometric analysis [29]. Briefly, the reagent was prepared by mixing 10 mmol 2,4,6-tris(2-pyridyl)-s-triazine (TPTZ)/L reagent with 20 mmol/L ferric chloride in acetate buffer (pH 3.6). Quantitative analyses were performed by the external standard method using ferrous sulphate (0.2 mmol/L) as the reference standard and correlating the absorbance (λ = 593 nm) with the concentration. 0.1 mL of beer/fruit juice sample was mixed in polystyrene cuvettes with 0.9 mL of distilled water and 3 mL of ferric complex. The data were expressed as Trolox equivalent (TE) of antioxidative capacity per 1 L of the beer (mmol TE/L). All measurements were performed in triplicate.

Free-Radical-Scavenging Ability by the Use of an ABTS +• Radical Cation
Another method used to measure the antioxidant activity of beers was the ABTS +• radical cation assay [30]. Sample of beer or fruit juice (0.03 mL) were mixed with 3 mL of ABTS +• solution with measured absorption of 0.700 at a wavelength of 734 nm. After 6 min the absorbance of samples was measured. Each sample was tested in triplicate. The data were expressed as mmol Trolox equivalent (TE) of antioxidative capacity per 1 L of the beer (mmol TE/L).

Adsorption of Volatile Compounds Using Solid Phase Microextraction (SPME)
Centrifuged and degassed beer (4 mL) was added to a 20 mL glass vial. 5 µL of 100 mg/100 mL aqueous emulsion internal standard solution (2-undecanone) was added. A magnetic stir bar was placed in the vial which was then closed with an aluminium membrane. The vial was placed on a IKA RCT basic magnetic stirrer (Guangzhou, China). The SPME holder needle with three-component universal fibre for SPME (DVB/CAR/PDMS 50/30 µm) was used to puncture the membrane. Heating to 40 • C and stirring was started and fibre was extended and held over the beer for 20 min. An analogical approach was used for each type of beer. The adsorption of volatile compounds contained in peeled hawthorn fruits added to beer was performed in order to find out which chemical compounds are characteristic for hawthorn fruit. Adsorption of compounds isolated from the hawthorn fruit was carried out in an analogous manner to the adsorption of volatile compounds from beer, but 5 g of hawthorn fruit instead of 4 mL of beer was added to the vial. The analysis of juice squeezed out of hawthorn fruit was carried out in an analogous manner to beer (4 mL juice volume). Measurements for each type of beer were performed in duplicate.

Gas Chromatography and Mass Spectrometry of Compounds Adsorbed on Fibre
Separation, identification and quantification of volatile compounds adsorbed on the fibre was carried out using a gas chromatograph connected to a Saturn 2000 MS Varian Chrompack mass detector (Palo Alto, CA, USA) with a ZB-5 column (Phenomenex, Torrance, CA, USA) (30 m length × 0.25 µm film thickness × 0.25 mm diameter).
Chromatographic conditions were carried out in accordance with the methodology of Calin-Sanchez et al. [31]. Scanning (1 scan/s) was carried out in the 35-400 m/z range using 70 mV electron ionization. The analyses were carried out using helium as a carrier gas with a flow rate of 1 mL/min using the following program for oven temperature: 40 • C at the beginning of the process, 5 • C/min up to 110 • C, 20 • C/min up to 270 • C. The initial temperature was maintained for 3 min. The injection port temperature of the chromatograph was 220 • C.

Sensory Analysis
The beers prepared in this study were subjected to an organoleptic assessment on a five-point scale using features such as clarity, aroma, colour, taste and overall impression. Beers were evaluated by a group of 15 people (21 to 32 years old), which consisted of 9 women and 6 men. Samples were given in coded plastic cups with a capacity of 200 mL. The temperature of the served beer was 8 • C.

Data Analysis
Volatile compounds separated from the beer were identified by comparing retention indexes (RI) with Kovats standards (KI exp. and KI lit.) and with NIST11 chemical standard libraries and by mass spectral analysis. Two standard matrices were also created (for the control sample chromatogram, i.e., beer without addition and hawthorn fruit). Chromatograms of the remaining samples (i.e., beer with addition of fruit and beer with addition of fruit juice) were integrated using retention time of compounds in those two standard matrices, using Mnova MS 12.0.1 software (Mestrelab Research, Santiago de Compostela, Spain). Data concerning content of volatile compounds in hawthorn juice and hawthorn fruit is available in Supplementary Materials (Table S1).
The results in this work were statistically analysed in the Statistica 12.5 program from Statsoft (Tulsa, OK, USA) using one-way ANOVA (α = 0.05). Differences between means were calculated using Duncan's test (α < 0.05).

Concentration of Carbohydrates, Glycerol, Ethanol and pH Value of Tested Beers
The ethanol content in the tested beers ranged from 35.43-41.63 g/L ( Table 1). The control beer without fruit addition (BC) had the highest content of ethanol, glycerol, maltose, maltotriose, as well as the highest pH value in relation to the other samples: beer with hawthorn fruit (BF) and beer with hawthorn juice (BJ). The addition of juice diluted and acidified the sample. This is indicated by a lower sugar content and lower pH characterising the samples with additives. BF compared to BJ, had a higher content of ethanol, glycerol, as well as maltose, maltotriose and dextrins.
The use of hawthorn fruit did not dilute the beer as much as juice addition. The fruit, after the primary fermentation process, was separated from the beer, so it did not significantly increase its final volume, in contrast to the samples with the addition of juice. This is reflected in the content of sugar, ethanol and glycerol per litre of beer, as well as in the pH value. Liu et al. analysed the sugar content in hawthorn fruits. The sugars identified in the fruits were: fructose, glucose, sucrose, sorbitol and myo-inositol [32]. The addition of hawthorn fruit or juice does not affect the content of maltose and maltotriose. The whole amount of these sugars was introduced into the beer with malt. The differences in sugar content between the samples are caused by the dilution of the samples by the addition of fruit or juice. Glucose content was not detected in any of the tested beers. Glucose is the sugar preferentially Foods 2020, 9, 775 6 of 14 used by the yeast Saccharomyces cerevisiae. The final product usually does not contain this sugar or contains a small amount of it [33]. The addition of fruit or fruit juices reduces the pH of the finished product. Nardini and Garaguso showed that the pH value for the classic ale beer style ranges from 4.39 to 4.73, while for fruit beers the pH reached even around 3.5 [18]. This was also confirmed in the studies by Kawa-Rygielska et al., where the effect of adding Cornelian cherry varieties on beer parameters was analysed. Control beer, without the addition, had a pH value of 4.59, while beers with the addition of juice had a pH in the range of 3.43-3.71 [20]. The content and composition of sugars and acids in the product is an important factor shaping the quality of fruit, and thus also fruit products, by affecting the taste and reception of the product by consumers [34].

Concentration of Total Polyphenols and Antioxidative Activity
The use of hawthorn fruit or juice resulted in a significant increase in the polyphenol content as well as the antioxidant activity tested by FRAP, ABTS +• and DPPH • methods. The addition of fruit caused an increase in the content of polyphenolic compounds in beer by 1.4 times (279.6 mg GAE/L), while the use of fruit juice more than twice (410.0 mg GAE/L) ( Table 2). A smaller increase in the polyphenol content was observed in BF compared to BJ. This can be explained by the poorer extraction of polyphenolic compounds from fruit into solution. Similar trends were observed for the beer's antioxidant activity. The addition of juice increases the total antioxidant capacity (ABTS +• ) 2.2 times, while the addition of fruit only 1.4 times. The ability to reduce the DPPH • radical solution increased the in BJ up to 6 times while in BF the increase was much less noticeable (only about 1.2 times). The FRAP assay showed a 2.6 times increase in antioxidant activity for BJ, and 1.7 times BF compared to BC. The addition of fruit or fruit juices in beer production is a method for significantly increasing the content of bioactive polyphenolic compounds in product. In the studies of Nardini and Garaguso, classic ale beers showed total polyphenol content in the range of 383-482 mg GAE/L, while beers with fruits (cherry, orange peel, grape, plum, raspberry, peach, apricot or apple) from 399 to up to 767 mg GAE/L [18]. The addition of fruits or fruit juices to beer can almost double the content of these biologically active compounds. BC contained 200.5 mg GAE/L polyphenolic compounds.
Fruit and vegetable peels can also be an interesting addition to the production of beer. Beers with eggplant skin extract were tested in a study by Horincar [22]. Polyphenol content in beers assessed in this study ranged from 0.443 to 0.610 mg GAE/mL, antioxidant capacity measured by DPPH • assay was 1.287-1.306 mmol TE/mL, while ABTS +• 0.095-0.107 mmol TE/mL. The content of polyphenols was of similar order as polyphenol content in BJ, but the antioxidant capacity of beer with addition of eggplant peel extract was significantly higher than antioxidant capacity of BJ. The main reason for such differences in tested beers was the amount of the addition used. The use of peels may be an interesting object for further research.
In addition to increasing antioxidant activity, polyphenolic compounds can affect the sensory characteristics of beer. They can introduce unpleasant bitter taste to beer and participate in the creation of turbidity [35]. Beer styles vary in polyphenol content. The type and dose of malts, hops, and the  [2] tested 12 commercially available beers including pilsner, lager, dark and black beers. The content of polyphenolic compounds ranged from 464.34 to 855.45 mg/L. Antioxidant activity of these beers was tested by DPPH • assay and ranged from 0.4 to 0.61 mmol TE/100 mL. These values were higher than those obtained by fruit beers in our experiment. In a study by Kawa-Rygielska et al. [20] on beers with the addition of Cornelian cherry juice, it was shown that Cornelian cherry juice increased the content of polyphenols in a greater extent than the addition of hawthorn juice. However, another study on non-alcoholic beers with the addition of Cornelian cherry juice, carried out by Adamenko et al. [36] showed a smaller increase in the antioxidant activity and the content of polyphenols in beers than the increase in the antioxidant activity and content of polyphenols in beer with hawthorn juice, but greater than determined in beer with fruits.

Concentration of Volatile Compounds in Tested Beers
The main volatile chemicals in alcoholic beverages are ethanol and carbon dioxide. Other compounds, such as esters, play a crucial role in creating beer aroma, despite being present at very low concentrations [37]. We found almost no research about volatile compounds in beers with fruit addition [38] and no works in which authors compared differences in volatile composition of beer made with addition of fruit or juice. In BC, 45 volatile compounds were identified ( Table 3). The largest groups among them were esters (20 compounds), sesquiterpenes (8 compounds), and alcohols (6 compounds). The greatest number (51) of volatile compounds was identified in BF. Like in the BC, the largest groups among the constituents were esters (25 compounds), sesquiterpenes (8 compounds) and alcohols (8 compounds). The analysis, in the case of BF, also revealed two peaks recognized as volatile constituents, which could not be identified (the mass spectra of unidentified constituents are available in the Supplementary Materials Figures S1 and S2). Other volatile compounds identified in BF, BJ and BC belonged to such groups as aldehydes, monoterpenes, ketones, organic acids and hydrocarbons. In a study by Gao et al. [39] about volatile compounds in hawthorn fruits, 61 volatile compounds were identified. Although more volatile constituents could be identified in fruit, our study managed to characterise more of the esters. The main reason for this phenomenon is probably yeast metabolism. Yeast cells can form esters by enzymatic chemical condensation of alcohols and organic acids [40]. The content of compounds identified in hawthorn fruit, such as ethyl butyrate, 1-hexanol, ethyl hexanoate, α-terpineol, ethyl octanoate, ethyl (4E)-4-decenoate and ethyl trans-2-decenoate was determined in all tested beers. The largest amount of compounds characteristic for hawthorn was recorded in BF and the smaller (except ethyl octanoate) in BJ. It is also worth noting that these compounds occur in a smaller amount, in the control sample. Therefore, they may also come from malt, hops or be products of yeast metabolism [41]. The addition of hawthorn fruits or hawthorn juice to beer introduced only six new compounds (acetic acid hexyl ester; propanoic acid 2-methyl, 2-methylbutyl ester; limonene, ethyl 5-methylhexanoate; octanoic acid methyl ester; isopentyl hexanoate). However, there were increases in the amount of almost all volatile compounds characteristic for the beer. The compound derived from hawthorn, which was found in the largest amount in BF, was α-terpineol, characterised by resinous and citrus aromas [42,43]. This constituent was found in the research about craft beers produced by mixed fermentation by Torulaspora delbrueckii non-conventional yeast and typical ale yeast [44]. In the control sample fermented solely with brewer's yeast, α-terpineol was not detected. In the samples fermented with Torulospora culture, the amount of α-terpineol was 3 times higher than in BJ, but 4.5 times lower than in BF. As in the study by Canonico [44], the amount of α-terpineol in BC was minuscule (0.1019 µg/100 mL). Another compound, one of the esters characteristic to the hawthorn fruit, was ethyl hexanoate. Among the esters characteristic of hawthorn fruit aroma, the concentration of ethyl hexanoate in BF was the highest. It is a compound with a sweet, fruity aroma [45,46]. Its presence was demonstrated in all tested samples: BF (11.2222 µg/100 mL), BJ (1.3959 µg/100 mL) and BC (0.9681 µg/100 mL). It is a compound which was characterised in a study about dotted hawthorn by Kucharska et al. [46]. In dotted hawthorn, ethyl hexanoate constituted 26.5% of all volatile components, but in beer, ethyl hexanoate constituted only 3.9% of all volatiles. In a study about wines made from hawthorn fruit juice, conducted by Zhang et al. [47], it was shown that concentration of ethyl hexanoate in wines was 7 times higher than in BJ, 10 times higher than in BC, and similar to the concentration in BF. Another ester identified in hawthorn fruit and all tested beers was ethyl (4E)-4-decenoate, which is characterised by fruit (pineapple), cognac and wax flavours [48,49]. It was detected in the tested samples in the amount of 0.1267 µg/100 mL for BJ, 0.4505 µg/100 mL for BF and 0.0985 µg/100 mL for BC. Concentration of ethyl octanoate was highest in BJ (8.3716 µg/100 mL) and lowest in BC (0.0656 µg/100 mL). This compound is characterised by a fruity aroma and can originate from hawthorn fruit as well as yeast metabolism. Ethyl octanoate has been found in various beers analysed by other authors, such as in beer with sorghum addition made by Chenge et al., where its determined concentration was two times higher than in BJ [41,46]. A compound not belonging to the group of esters, which was identified in the fruits of hawthorn, was 1-hexanol, which is a compound from chemical family of alcohols. It has a sweet, alcoholic aroma [50,51]. Its highest concentration was found in BF (0.0819 µg/100 mL). BJ contained 0.0327 µg/100 mL of this compound, while BC contained 0.0164 µg/100 mL. This compound was also found in greater concentration in fruit beers with addition of quince than in control beer in study of Zapata et al. [38]. The most abundant sesquiterpenes in beers were humulene and humulene epoxide I. Concentration of humulene was the highest in BF (1.9659 µg/100 mL), and the lowest in BJ (0.4772 µg/100 mL). Meanwhile, humulene epoxide I was found in BF in the amount of 3.1127 µg/100 mL and BJ in the amount of 0.2812 µg/100 mL. In BC, the concentration of these compounds was greater than in BJ (respectively 0.5118 µg/100 mL and 0.2969 µg/100 mL). In classic beers, these compounds come from hops [52], while this study may indicate that these sesquiterpenes can be found in hawthorn fruit, but only in the flesh (because in BJ their quantity is smaller than in BC). Another compound, β-farnesene, with the lowest concentration in BJ, which is naturally found in classic beers, also comes from the sesquiterpene group [53]. It is a compound characterised by woody and citrus aromas [54,55]. Its concentration in BJ was 0.0136 µg/100 mL, while in BC it was 0.0492 µg/100 mL, and in BF 1.6383 µg/100 mL. The total concentration of all volatile compounds in the tested beers was the highest for BF which was 287.8529 µg/100 mL, which was 6.5 times higher than in BJ, which contained 43.8884 µg of volatile compounds per 100 mL of beer, and 8 times higher than BC, which contained 35.9348 µg of volatile compounds per 100 mL of beer.

Organoleptic Analysis
The aroma, colour, clarity, foaminess and overall impression of beer were assessed. The beer that achieved the highest rating in terms of aroma, taste, clarity and overall impression was BF (Table 4, Figure 1). This is in line with the results of research carried out by Adamenko et al. [36] on beers with the addition of Cornelian cherry, where beers with the addition of fruit also received better marks in terms of taste than control beers. Statistical analysis showed no significant differences between BC and BF clarity, although the addition of juice negatively affected this attribute. The use of whole fruit caused a significant deterioration in beer foaminess. Control beer got the highest grade for foaminess; however, statistical analysis did not show a significant difference between the BC and BJ results. The addition of juice does not worsen the beer foam quality. It is possible that a high concentration of polyphenols in beer can have a detrimental effect on beer foaming [2]. Although some authors show a positive effect of high polyphenols on beer foam stability [56]. The results obtained are contradictory. It has not been clearly defined how polyphenolic compounds affect beer foam. BF achieved the best note for taste criterion. In this aspect of sensory analysis, it obtained significantly higher results than other analysed beers. BJ and BF beers received better aroma ratings than BC. This may be due to a higher concentration of esters characterised by floral and fruity notes (ethyl hexanoate or ethyl octanoate). Zapata et al. [38] presented similar results. They analysed the composition of the quince beer aroma. The study showed that the higher concentration of ethyl esters resulted in better scores in sensory analysis.   significantly higher results than other analysed beers. BJ and BF beers received better aroma ratings than BC. This may be due to a higher concentration of esters characterised by floral and fruity notes (ethyl hexanoate or ethyl octanoate). Zapata et al. [38] presented similar results. They analysed the composition of the quince beer aroma. The study showed that the higher concentration of ethyl esters resulted in better scores in sensory analysis.

Conclusions
The study indicates that the addition of hawthorn (Crataegus punctata) makes it possible to obtain beer with an increased total content of polyphenolic compounds and antioxidant activity. Beer without hawthorn addition had the highest content of ethanol carbohydrates, glycerol and highest pH value. The total content of polyphenolic compounds was affected by both the addition of fruit and the form of the fruit added. Hawthorn in the form of juice made it possible to obtain a much higher content of bioactive polyphenolic compounds in beer than addition of fruit. Nevertheless, beer with fruit addition was characterised by the greatest amount of volatile compounds, eight times higher than in the control sample. Enrichment of the beer by hawthorn fruit mostly increased the concentration of the volatile compounds characterised by fruity and sweet aromas. The results of sensory analysis indicate that addition of hawthorn fruit results in an improvement of such characteristics as taste, aroma, clarity and overall impression in a higher degree than hawthorn juice addition. Hawthorn fruit and its juice can be used as complementary raw material in the production of beer to increase its biological activity and improve its taste and aroma. It can also contribute to greater consumer interest in the product.