Broomrape as a Major Constraint for Grass Pea ( Lathyrus sativus ) Production in Mediterranean Rain-Fed Environments

: Grass pea ( Lathyrus sativus ) is an annual legume crop that is currently underutilized but has the potential for reintroduction into Mediterranean rain-fed farming systems. In this study, we compared the adaptation of breeding lines in multi-environment ﬁeld testing, which had wide variation for precocity, grain yield and broomrape infection. Heritability-adjusted genotype plus genotype-by-environment interaction (HA-GGE) biplot and non-metric multidimensional scaling (NMDS) were performed to determine the e ﬀ ect on genotype (G), environment (E) and G × E interaction on grain yield, precocity and broomrape infection. Precocity was associated with reduced broomrape infection, and this with increased grain yield. Step-wise regression analysis revealed that the broomrape infection had the highest inﬂuence on grain yield, whereas precocity had a lower e ﬀ ect. Rain and humidity and mild temperatures before and during ﬂowering were the climatic factors most inﬂuential on broomrape. Accessions with a shorter growth cycle su ﬀ ered lower broomrape infection and were more productive in the environments with a high broomrape incidence. Accessions with longer growth cycle su ﬀ ered overall higher broomrape infection and were therefore more productive in the environments with low or moderate broomrape incidence.


Introduction
Grass pea (Lathyrus sativus L.) is an annual temperate legume widely grown in South Asia and sub-Saharan Africa, but which is in regression in Mediterranean environments [1]. Common experience in Mediterranean countries is a severe decrease in grass pea cultivation during the last century, with little or negligible trade. For instance, in Spain, the acreage decreased from 250,000 ha in 1925 to 636 ha in 2018 [2,3]. A major factor for this was the ban for human consumption in some countries due to the occurrence of lathyrism, a neurodegenerative disease in humans and domestic animals [4,5]. However, the ban is being reassessed due to the current understanding that the neurotoxin β-ODAP content is highly influenced by climatic and edaphic conditions [6][7][8], and that lathyrism disease is caused only when there is an overconsumption of grass pea in a non-balanced diet, that can be prevented by the addition of sufficient cereals, fruits and vegetables [9]. Therefore, for instance, this ban has been removed recently in Spain [2]. Variation for β-ODAP content is available in grass pea which is considered to be lower in Mediterranean white seed types, although also selected into colored seeds [6,10,11].
Reducing β-ODAP would be desirable, however, there is no particular risk in the occasional intake of grass pea as part of a balanced diet [2,4]. In this study, we focused on the adaptation and grain yield, which should be improved in order to make the crop's cultivation feasible. Due to its rusticity and potential in marginal environments, there is a renewed interest to re-introduce the grass pea into Mediterranean rain-fed cropping systems. In Mediterranean regions, the root parasitic weed crenate broomrape (Orobanche crenata Forsk.) is known as one of the major constraints for grass pea cultivation [12,13]. Therefore, efforts have been made to gather and characterize the adaptation of landraces and to submit them to breeding in order to exploit the potential of the species [4,6,12,[14][15][16][17].
The present study aimed at evaluating the yield performance and stability among grass pea accessions in Mediterranean rain-fed cropping systems.

Plant Materials and Experimental Design
The performance of 8 grass pea breeding lines was studied at 17 location-year environments ( Table 1). These lines were selected from previous field studies [12,13], gathering among the best yielders in those studies, four lines with colored flowers and seeds (Ilat3, Ilat10, Ilat11 and Ilat18) and four with white and large seeds (Tblanco, Tpinto, Titana and Lasana), which are preferred by Mediterranean consumers [11]. At each location, a randomized complete block design with three replications was used. The experimental unit consisted of 3 m-long rows per accession separated by 0.35 m, 10 plants per row. Seeding was performed manually by mid-December at each site-year. Weeds were controlled by hand weeding, with no herbicides or pesticides applied. Days to flowering (DtF) was estimated in five environments by weekly recording the date in which 50% of the plants of each plot had at least one fully opened flower. Numbers of emerged broomrape plants per row were recorded and referred to as the number of broomrapes per grass pea plant. The presence of naturally occurring pest and disease was also recorded, estimating the percentage of canopy coverage or seed damage. The plots were harvested manually at full maturity by late April, early May, depending on the environment. Seeds were then threshed and the grain yield assessed.

Statistical Analysis
A combined ANOVA for randomized complete-block designs within each year-location environment was carried out using SAS ® 9.3 (SAS Institute Inc., Cary, NC, USA). F ratios, used to test effects for randomized complete block experiments combining location-year environments (Table 1), were determined according to McIntosh [18]. Prior to each ANOVA, tests for the normality and equality of variance were conducted for each dependent variable.
To determine the stability and identify superior accessions across environments, heritability-adjusted genotype plus genotype × environment interaction (HA-GGE) biplot analyses were conducted [19,20], comprising the seasons from 2008 to 2018, since it takes into consideration any heterogeneity among environments by giving weights to the test environments proportional to their root square heritability. As the environments have different heritabilities (data not shown) for the same trait, the HA-GGE biplot was most appropriate for the visual evaluation of the test environments and genotypes [19]. Analyses were made with the SAS ® 9.3 (SAS Institute Inc.) program developed by Burgueño et al. [21] to graph GGE biplots.
The genotype × environment (G × E) two-way tables were first centered with the respective means for the environments, multiplied by √ H and then divided by the SD of the respective environment [19].
The HA-GGE biplot shows the first two principal components (PC1 and PC2) derived from the previous two-way table of yield to the singular value decomposition [22,23].
Singular value partitioning is achieved by providing a scaling factor f to obtain alternative accessions and environment scores. We chose the most straightforward variant called symmetric scaling (f = 0.5) since it bears most of the properties associated with other scaling methods [24].
The target environment axis abscissa (TEAa) is represented by a straight line drawn through the biplot origin and the average environment, which was defined by the mean ordinates of all environments in the biplot. The main effects of genotypes (G) are represented by the projection of genotypes onto this axis. These projections provide us with the contribution of each genotype to G, so genotypes may be ranked along the TEAa, with the arrow pointing to a higher mean yield. The TEA ordinate (TEAo) shows the contribution of each genotype to the interaction G × E, thus giving information about the genotypic stability or instability (consistency or inconsistency across environments). The best genotype would be that with the highest yield (higher positive projection on TEAa) and the highest stability, i.e., projection on TEAo close to 0 [22]. Similarly, the ideal environment would be the one showing a high projection value onto the TEA abscissa (further discriminating the principal effects of genotypes) and a small absolute projection value onto TEAo (more representative of all the tested environments) [22].
To evaluate the influence of environmental factors on the broomrape infection, 21 climate variables (Supplementary Table S1) were subjected to non-metric multi-dimensional scaling ordination (NMDS) [25]. These climate variables were obtained from the Junta de Andalucía (https://www.juntadeandalucia.es/agriculturaypesca/ifapa/ria/servlet/FrontController) and included the maximum, minimum and average temperature, maximum, minimum and average humidity, accumulated radiation, evapotranspiration and accumulated rain during pre-flowering, at flowering and post-flowering period. To decrease the probability that the result of the NMDS analysis would reflect a local stress minimum rather than the overall minimum, we repeated the NMDS analysis 20 times, each time starting from a different random configuration, and selected the two-dimensional solution with the lowest stress. Analysis was made using the PAST software [26].
Stepwise linear regression was applied to describe the impact of broomrape, of flowering date, and all the climatic parameters on yield. Significance levels for introducing and deleting a variable was set at α = 0.05.

Analysis of Variance
The combined analysis of variance revealed that all main effects (environments (E), genotypes (G) and G × E interaction) were statistically significant for the grain yield, number of broomrapes per plant and flowering date ( Table 2). The environment explained 52 to 81% of the total variation (G + E + G × E sum of squares), G and G × E accounted 3 for 29% and 11 to 18%, respectively. Grain yields of the accessions for each environment are shown in Table 3 and Figure 1. The average yield per environment over accessions was the highest at TOM09 (4675 kg ha −1 ) and ALM17 (4550 kg ha −1 ) and lowest at ESC10 (209 kg ha −1 ). The average grain yield over environments and accessions was 2143 kg ha −1 for L. sativus, the highest for Titana (2643 kg ha −1 ) and Lasana (2517 kg ha −1 ) and the lowest was for Tblanco (1714 kg ha −1 ) and Tpinto (1865 kg ha −1 ).   Table 1. TEAa is the target environment axis abscissa and TEAo is the target environment axis ordinate.

HA-GGE Analysis
The first two PCs of the HA-GGE model for yield, number of broomrapes per plant and days to flowering explained 77, 93 and 82%, respectively, of total G + G × E. On the other hand, (G + G × E)/(E + G + G × E) yielded a values of 19, 30 and 47%, respectively (Figures 1-3). This fulfilled the requirements of Yang et al. [24], who established that for a biplot to be useful, the first two PCs should be higher than 60% and the (G + G × E)/(E + G + G × E) ratio should be higher than 10%.  Table 1. TEAa is the target environment axis abscissa and TEAo is the target environment axis ordinate.

HA-GGE Analysis
The first two PCs of the HA-GGE model for yield, number of broomrapes per plant and days to flowering explained 77%, 93% and 82%, respectively, of total G + G × E. On the other hand, (G + G × E)/(E + G + G × E) yielded a values of 19%, 30% and 47%, respectively (Figures 1-3). This fulfilled the requirements of Yang et al. [24], who established that for a biplot to be useful, the first two PCs should be higher than 60% and the (G + G × E)/(E + G + G × E) ratio should be higher than 10%.  Table 1. Environments ALM17, ALM18, ESP17, ESP18 were free of broomrape infection and therefore are not included in this biplot.  Table 1. Environments ALM17, ALM18, ESP17, ESP18 were free of broomrape infection and therefore are not included in this biplot.
Agronomy 2020, 10, x FOR PEER REVIEW 6 of 12  Table 1. Environments ALM17, ALM18, ESP17, ESP18 were free of broomrape infection and therefore are not included in this biplot. The yielding ability and the stability of an accession is shown in the biplot (Figure 1) by the "TEA". The projection of an accession over the average environment axe (TEAa) indicates its mean performance across all environments, and its projection over the TEAo indicates its stability [27]. Figure 1 shows that Titana, and Lasana were the highest yielding accessions on average (high absolute primary scores, TEAa). Tblanco and Tpinto were the lowest yielding and most stable over the environments (low absolute secondary scores, TEAo) [27].
Among biotic constraints, our study showed broomrape as the most significant one. Little incidence of any other pest or disease was observed, with the exception of some powdery mildew appearing late in the season, with disease severity lower than 20%, and bruchus weevil seed infestation lower than 20% (data not shown). These values are lower than those reported in pea in those locations (up to 100% for powdery mildew and 80% for weevil) [28]. The global average of broomrape infection over accessions and environments was 0.69 broomrapes per grass pea plant (Table 4). Average infection over environments was lower for the Ilat accessions (range of 0.35-0.47 broomrapes/plat) and highest for Titana and Tblanco (>1 broomrapes/plant).  Ilat3, Ilat10, Ilat11 and Ilat18 are more distantly located on the left to TEAo vector (lowest infection), whereas Titana and T blanco are more to the right (highest infection) ( Figure 2). The response of the four Ilat accessions was rather stable among environments, as shown by their shorter projection on TEAo (the closer to the TEAa vector). On the contrary, the response of Lasana and Tblanco was less stable, with Lasana being more infected in some environments and Tblanco in others.
The earlier accessions (Table 5) were the Ilat accessions (88 to 98 DtF), with the remaining accessions having >104 DtF. This is shown in the biplot (Figure 3) with the earlier accessions to the left of TEAo, so the DtF increase as we move to the right. The closer they are to TEAa, the more stable the trait, which was more affected by environment for Ilat3 and Lasana.

Effect of Broomrape Infection, Flowering Date and All Climate Variables on Yield
The produced regression model is based on 21 climatic parameters, with the broomrape infection and flowering date as the independent variables. The dependent feature refers to the yield of grass pea (kg ha −1 ). The most appropriate model obtained in the step-wise regression analysis was:

Effect of Broomrape Infection, Flowering Date and All Climate Variables on Yield
The produced regression model is based on 21 climatic parameters, with the broomrape infection and flowering date as the independent variables. The dependent feature refers to the yield of grass pea (kg ha −1 ). The most appropriate model obtained in the step-wise regression analysis was: This result is in agreement with the coefficients of correlation, so the coefficient above commented between the broomrape infection and grain yield (r = −0.78; p < 0.0001) and between maximum humidity at flowering and grain yield (r = −0.48; p < 0.0001).

Discussion
Great morphological variation is reported in grass pea, showing a clear grouping in two major types, one with a covering blue flower and colored seeds that is typical of the Indian subcontinent, and other types with white flowers and seeds that have a more western distribution [29]. Ilat studied accessions, originating from ICARDA-CGIAR (International Center for Agricultural Research in the Dry Areas) programs, that happened to fit into to the first group, having blue flowers and colored seeds, whereas Tblanco, Tpinto, and Titana, deriving by selection among Iberian landraces, and Lasana, deriving from a landrace of unknown origin, fitted the second group, with mostly white flowers and larger and whiter seeds. This is the seed type preferred by Mediterranean consumers [11,29]. We found that these accessions have a later phenology, which is in line with previous reports [11,12].
Grass pea breeding efforts in past decades focused on reducing the β-ODAP content. However, the long-term results of these efforts are often questioned as the β-ODAP content is highly influenced by climatic and edaphic conditions and lathyrism is caused only when there is an overconsumption of grass pea in a non-balanced diet (2,(4)(5)(6)(7)(8), which made breeders reconsider the yield stability as a higher concern. We therefore focused on the adaptation and grain yield, which should be improved in order to make the crop feasible. We found wide variation for precocity, grain yield and broomrape infection, which have major implications in the crop performance.
In Mediterranean regions, the broomrape is acknowledged as a major constraint for grass pea cultivation [12,13]. Broomrape is also a constraint to other grain and forage legumes [30], but the grass pea seems to be particularly sensitive, suffering a higher yield penalty than other legumes [31]. Some levels of resistance have been reported in grass pea [13,32], which might be the result of a combination of different escape and resistance mechanisms [33]. We observed an association between late phenology and broomrape infection, which is in agreement with previous reports of grass pea [12,13] and other legumes [34][35][36][37]. In addition to the host, the infection severity of broomrape strongly depends on parasitic seedbank density and on environmental factors such as temperature and rain. We found that rain and mild temperatures before and during flowering were the climatic factors most influential on broomrape infection. These effects have been well documented, with cool winters and limited rain reported to reduce infection [34,38,39].
In agreement with previous reports [12], we found broomrape infection to be reduced with grass pea precocity, and this with increased grain yield. Accessions with a shorter growth cycle suffered lower broomrape infection and were more productive in the environments with high infection. Those with a longer growth cycle (Titana, Lasana, Tpinto and Tblanco) suffered overall a higher broomrape infection and were therefore more productive in the environments with low or moderate broomrape infection. However, step-wise regression revealed that is reduced broomrape infection which has the highest influence on grain yield, with precocity alone having very little effect. This might be due to a compensation of the beneficial effect of grass pea precocity facilitating the escape from broomrape and drought [13], with the detrimental effects of reducing the crop cycle and limiting the potential yield [7,16,17].
We concluded that the grass pea has potential for reintroduction into Mediterranean rain-fed farming systems. Broomrape infection appears as a serious constraint, having a major influence on grain yield. Early cultivars are recommended for areas prone to high broomrape infection. Cultivars with a longer growth cycle suffer higher broomrape infection but can be more productive in the environments with low or moderate broomrape incidence.