Genetic and Genomic Diversity in a Tarwi ( Lupinus mutabilis Sweet) Germplasm Collection and Adaptability to Mediterranean Climate Conditions

: Lupinus mutabilis (tarwi) is a species of Andean origin with high protein and oil content and regarded as a potential crop in Europe. The success in the introduction of this crop depends in part on in depth knowledge of the intra-speciﬁc genetic variability of the collections, enabling the establishment of breeding and conservation programs. In this study, we used morphological traits, Inter-Simple Sequence Repeat markers and genome size to assess genetic and genomic diversity of 23 tarwi accessions under Mediterranean conditions. Phenotypic analyses and yield component studies point out accession LM268 as that achieving the highest seed production, producing large seeds and e ﬃ ciently using primary branches as an important component of total yield, similar to the L. albus cultivars used as controls. By contrast, accession JKI-L295 presents high yield concentrated on the main stem, suggesting a semi-determinate development pattern. Genetic and genomic analyses revealed important levels of diversity, however not relatable to phenotypic diversity, reﬂecting the recent domestication of this crop. This is the ﬁrst study of genome size diversity within L. mutabilis , revealing an average size of 2.05 pg / 2C (2001 Mbp) with 9.2% variation (1897–2003 Mbp), prompting further studies for the exploitation of this diversity.


Introduction
The genus Lupinus includes more than 280 species [1], approximately 90% of which are native and widely distributed throughout the American continent [2,3], with greater inter-and intra-specific genetic variability than in Euro-African species. Lupinus mutabilis Sweet (also known as tarwi, chocho, altramuz and Andean lupin) is native from the Andean region in South America. The species is autoand allogamic with wide variability of flower, stem and seed colours, and exhibits indeterminate growth [4,5]. It has been domesticated in the Andean region and used for grain production, forage, green manure, fixing atmospheric nitrogen and soil conservation [4,6]. In spite of their high alkaloid content [7,8], tarwi seeds have high nutritional value, containing up to 53% protein and 24% lipids [9]. The nutritional attributes found in tarwi are supposedly better than those in soybeans [10] and for this reason it is called Andean soybean [11]. Tarwi protein is rich in globulins (43%-45%) and albumins (8%-9%) and the oil has high quality and does not require industrial removal of the linolenic acid like in soybean [12,13]. Additionally, low alkaloid (<0.1%) lines have been selected in L. mutabilis [14].
Tarwi exhibits key traits of domestication, including indehiscent pods and seeds with permeable tegument, representing a locally important crop in several Andean areas [15]. Recently the species Agronomy 2020, 10, 21 3 of 24

Plant Materials
A total of 23 L. mutabilis accessions were selected from the ISA Lupinus germplasm including five accessions provided by the Julius Kühn-Institut (JKI), Germany (e.g., Table 2). Lupinus albus cultivars Misak and Mihai were used as reference in the morphological characterization because of their high adaptation to the Mediterranean conditions and as outgroups/standards in the ISSR marker and genome size analyses.

Morphological Analysis
Field experiments were conducted at Tapada de Ajuda in Lisbon (coordin: 38.709133, −9.182976, alt: 60 m) on a vertisol in the 2016/17 (sowing date: 29 December) and 2017/18 (sowing date: 18 December) seasons under rain-fed conditions. Meteorological data were collected daily from the weather station located adjacent to the field. Soil water balances were calculated according to Allen et al. [39].
The experimental design adopted was randomized block with three replicates. Each replicate was composed of 26 1.8 m 2 -plots with 20 plants in each plot (immediately surrounded by a 60 cm-wide edge of L. albus 'Misak' plants to avoid border effects) and the total number of plots in the assay was 78. For morphological characterization, 10 plants of each plot were selected as recommended by Talhinhas et al. [40].
Data of morphological characterization were obtained based on Lupinus spp. descriptors [41], as listed in Table 1. Yield components and vegetative traits were analysed considering a two-factor experimental design (genotype and year), with differences being statistically analysed using the Kruskal Wallis test. Characteristics for multivariate analysis were selected based on correlation coefficients and heritability values [40]. Variables with correlation above 0.85 were considered redundant and thus one of them was excluded. Meanwhile, variables with low heritability (<65%) were also excluded, as these were explained by environmental factors.
Univariate analysis (UA) was performed to compare each individual characteristic across the accessions. Before running the UA, normality and homogeneity of variances was tested. Since data did not follow normal distribution and the variance was not homogeneous, an analysis of variance (ANOVA) based on rank transformation for non-parametric analysis was performed [42]. Post-hoc Tukey's honest significant difference (HSD) test of means was performed for all variables at 5% significance. Afterwards, broad sense heritability (H 2 ), genotypic variance (σ 2 g ), phenotypic variance (σ 2 P ), phenotypic coefficient of variance (PCV) and genotypic coefficient of variance (GCV) were estimated to understand the genetic variation between accessions and environment, as well as the genetic effects on different traits, following Mazid et al. [43].
Multivariate analysis was performed for all 25 accessions and all characteristics selected and represented in a single graphic, as described by Talhinhas et al. [40]. Standardized morphological data transformation (mean = 0, and standard deviation = 1) was performed before conducting multivariate analysis. Cluster analysis was performed based on Euclidean distance and average method for the 25 accessions. A dendrogram was constructed using an unweighted pair group method of arithmetic mean (UPGMA) algorithm. Principal component analysis (PCA) was performed and eigenvectors and eigenvalues were projected to visualize the components. All analyses were performed in the RStudio program version 1.1.456 (The R consortium, Boston, USA). Total thousand seeds weight = TW/TNS × 1000 (g) SWBLR Seed weight/total branch length ratio = TW/TBL × 100 (g/m) 1 Characteristics related with secondary and tertiary branches were excluded due to insufficiency of data; 2 Redundant or non-independent characteristics excluded of multivariate analysis based on the correlation coefficient (r > 0.85); 3 Characteristics excluded of multivariate analysis due to presenting low value of heritability (<0.65); 4 number; 5 Metric measurement; 6 dim.-dimensionless.

Molecular Analysis
Young but fully expanded leaves of the 23 L. mutabilis accessions and of the two L. albus reference cultivars were collected and immediately frozen in liquid nitrogen and stored at −80 • C. Freeze-dried vegetal material was used for DNA extraction using the DNeasy ® Plant mini kit (Qiagen, Hilden, Germany) according to the manufacturer instructions. The DNA quality and quantity were estimated using spectrophotometry in the Gen5 program, and electrophoresis using a 1% agarose gel. The stock solution of DNA was diluted with sterilized water to make a working solution with a concentration of 10 ng/µL to be used in amplifications.
For molecular characterization, six ISSR primers were selected (Table 6) from those reported by Talhinhas et al. [44] based on preliminary analyses of a limited set of accessions. The Polymerase Chain Reaction (PCR) amplification for all primers was carried out under the following conditions: pre-denaturation 4 min at 94 • C, 40 cycles of 30 s at 94 • C, 45 s at 52 • C and 2 min at 72 • C, and a final extension at 72 • C for 10 min. The PCR reactions were performed in a final volume of 10 µL containing 20 ng of DNA, 0.5 µM of primer and 5 µL of dNTP + Taq DNA polymerase (NZYTaq II DNA polymerase, NZYTech, Lisbon, Portugal). After amplification, products were separated by electrophoresis in a 2% agarose gel stained using GreenSafe Premium (NZYTech).
The ISSR bands were scored in a binary matrix as presence (1) or absence (0) for each accession and for each fragment size. Based on the binary matrix, parameters such as percentage of polymorphic and monomorphic bands were determined and discriminatory power of primers was calculated based on the polymorphic information content (PIC), effective multiplex ratio (EMR), resolving power (RP) and marker index (MI). PIC value is the probability for detecting polymorphism by a primer or primers combination between two randomly drawn genotypes and can be calculated using the formula PIC = 1 − Σpi 2 , where pi is the frequency of occurrence of polymorphic bands in different primers [45]. The effective multiplex ratio was calculated using the formula EMR = npβ; where ß is the fraction of polymorphic markers and is estimated after considering the number of polymorphic loci (np) and non-polymorphic loci (nnp) as ß = np/(np + nnp) [46]. Marker index (MI) is the primer capacity to detect polymorphic loci among different genotypes and was calculated as EMRxPIC. Resolving power (RP) is the ability of primers to distinguish between genotypes and was calculated as RP = ΣIb, where Ib is the informative fragments and can take values of: 1 − [2|0.5 − p|]; p is the proportion of total genotypes containing the band [47]. Genetic similarity was obtained according to the Jaccard similarity index. The results were used for the construction of ISSR and morphological traits dendrograms, in order to evaluate the similarity relations between the genotypes. Dendrograms were constructed based on UPGMA grouping and the ISSR results were correlated with morphological traits.

Flow Cytometry
For each accession, young leaves in healthy conditions were randomly collected and immediately analysed in the laboratory. Nuclear DNA content was measured by flow cytometry. Solanum lycopersicum 'Stupické' (2C = 1.96 pg; [48]) was tested as DNA standard but its genome size showed to be too close to that of L. mutabilis. Therefore, we tested L. albus as DNA standard (2C = 1.20 pg; [49]) and for such L. albus 'Misak' was validated as standard by comparison to S. lycopersicum 'Stupické' and Raphanus sativus 'Saxa' (2C = 1.11 pg; [48]). Each L. mutabilis accession, together with the standard, was chopped with a razor blade in the presence of 1 mL of buffer (Woody Plant Buffer; [50]). The nuclear suspension obtained was then separated from plant debris using a 30 µm nylon filter. After filtration, 50 µg/mL of propidium iodide (PI; Sigma-Aldrich) were added to stain DNA and 50 µg/mL of RNase (Sigma-Aldrich) ware added to prevent staining of double stranded RNA. The samples were maintaining at room temperature and analyzed using a CyFlow Space flow cytometer (Sysmex, Norderstedt, Germany) equipped with a 30 mW green solid-state laser emitting at 532 nm for optimal PI excitation. The reproducibility of results were assessed using five independent replicates for each accession. FloMax software v2.4d (Sysmex) was used to measure nuclear DNA content and three graphics were generated from data measurement: fluorescence pulse integral in linear scale (FL); fluorescence pulse integral in linear scale versus time; and fluorescence pulse integral in linear scale versus side light scatter in logarithmic scale (SSC). The absolute DNA amount of a sample was calculated based on the values of the G1 peak means, as suggested by Doležel and Bartoš [51]: The results generated from 2C DNA (in picogram) were transformed to million base pairs using the following conversion: 1 pg = 978 Mbp [52]. Coefficient of variation (CV, %) of G1 peaks in the FL histograms, and estimates of the CV of the genome size of each accession were used to assess the reliability of the results. Intra-specific genome size comparison was carried out using Kruskal Wallis test (α = 0.05) because genome size data did not exhibit normal distribution. Data analysis was done in RStudio Program Version 1.1.456.

Morphological Characterization and Genetic Parameters among Accessions
Studies on the genetic variability are important because they generate relevant data for breeding programmes and can be used as basis for development and selection of superior genotypes. Here, we used morphological characterization and genetic parameters to evaluate the variability of a L. mutabilis germplasm collection under Mediterranean conditions.
Meteorological conditions during the trial were typical of the Mediterranean climate ( Figure S1), although rainfall was well below average during autumn and winter and above average during spring in 2017/18, while rain was scarce in April 2016. Two-way ANOVA based on rank transformation was performed and revealed that all morphological traits exhibited significant difference at p-value < 0.05. Tables 2 and 3 show the mean values, homogeneous groups and p value obtained in a two-factor experimental design for morphological and reproductive characteristics, respectively. Differences analysis results for morphological traits of each year are given in Supplementary Material (Table S1). The statistical analysis of the results depict those genotypes showing differences that were consistent over the two years.  The total seed weight (TW) per plant (Table 3) varied 3.7×, ranging between 3.7 g per plant (accession JKI-L210) and 13.8 g per plant (accession LM268), the latter attaining a projected productivity estimated at 1533 kg/ha, although less than half of the total yield of the L. albus accessions. Dissecting yield components evidences additional variability among the accessions ( Table  3). The total number of pods (TNP) per plant varied nearly 2.1×, with a maximum of 25 pods per plant for accession Potosi-ISA. The total number of seeds (TNS) reached a maximum of 67.9 seeds per plant (accession LM34). The average number of seeds per pod (TNSPP) is 2.7, ranging between 2.2 (accession LM268) and 3.3 (accession JKI-L210). The total thousand seeds weight (TTSW) attained a global average of 187.0 g, varying 2.85× between 101.4 g (accession JKI-L210) and 289.2 g (accession LM268). LM268 was the only L. mutabilis accession producing more yield on the primary branches than on the main stem (40% and 60% of total yield on the main stem and on primary branches, respectively), following a similar pattern to that of L. albus cultivars. The accessions JKI-L295 and JKI-L210 produced about 80% of seed weight on the main stem. Unlike accession LM268, several L. mutabilis accessions reached superior seed yields (over 10 g per plant) while concentrating over 60% The results presented in Table 2  The JKI-L309 accession presented low values for TBL, suggesting that this may be a semi-determinate genotype. Stem and flower colors varied among accessions, with no clear correlation to morphologic traits. Figure 1 depicts the four groups of flower and stem colors. 0.000 0.000 0.000 0.000 0.000 0.000 0.000 1 Homogeneous groups-accessions sharing the same letter for each trait are not statistically different; 2 SFC-stem and flower colour, according to Figure 1; 3 Full name of acronyms and description of the respective morphologic traits are given in Table 1; 4 p-value taking into account the accessions; 5 p-value taking into account the experiments; 6 p-value taking into account the interaction between accessions and experiments.
The total seed weight (TW) per plant (Table 3) varied 3.7×, ranging between 3.7 g per plant (accession JKI-L210) and 13.8 g per plant (accession LM268), the latter attaining a projected productivity estimated at 1533 kg/ha, although less than half of the total yield of the L. albus accessions. Dissecting yield components evidences additional variability among the accessions ( Table 3). The total number of pods (TNP) per plant varied nearly 2.1×, with a maximum of 25 pods per plant for accession Potosi-ISA. The total number of seeds (TNS) reached a maximum of 67.9 seeds per plant (accession LM34). The average number of seeds per pod (TNSPP) is 2.7, ranging between 2.2 (accession LM268) and 3.3 (accession JKI-L210). The total thousand seeds weight (TTSW) attained a global average of 187.0 g, varying 2.85× between 101.4 g (accession JKI-L210) and 289.2 g (accession LM268). LM268 was the only L. mutabilis accession producing more yield on the primary branches than on the main stem (40% and 60% of total yield on the main stem and on primary branches, respectively), following a similar pattern to that of L. albus cultivars. The accessions JKI-L295 and JKI-L210 produced about 80% of seed weight on the main stem. Unlike accession LM268, several L. mutabilis accessions reached superior seed yields (over 10 g per plant) while concentrating over 60% of their yield on the main stem: CM157, I82 and LM27. For the comparison of seed yield and vegetative development (Table 3), the seed weight/total branch length ratio (SWBLR) was calculated. SWBLR average was 2.0 g of seeds per meter of branch length in L. mutabilis (23.1 g/m in L. albus), ranging between 1.1 g/m (accession LM32) and 3.3 g/m (accession Mutal). In Table 4 are presented the average, phenotypic and genotypic variance with their respective coefficient of variation and heritability for the 2016/17 and 2017/18 experiments. Higher values of phenotypic and genotypic variances were observed for TSWPB, TSWMS, TBL, TNS, SLPB, NSPB and WSPPB (see Table 1   Correlation is an important test and is used to assess relationship and associations between variables and is frequently applied in several studies. The results generated after correlation coefficient analysis using the Spearman method for morphological traits in experiments of 2016/17 and 2017/18 are presented in Tables S2-1  From the correlation data, heat maps were constructed (Figures 2 and 3) using euclidean distances and the UPGMA method, where in the vertical columns are the clusters of morphological traits while in the horizontal lines are the clusters of accessions. Dark red colors represent lower values while the dark blue are higher values. Figure 2 corresponds to the heat map obtained from 2016/17 data and Figure 3 was obtained with 2017/18 data. Six groups of morphological traits could be drawn in Figure 2 and five in Figure 3. In both figures, group 1 is related to the reproductive capacity of pods, defined in the Figure 2 by the characteristics NSMS, SWBLR and TW and in Figure 3 by PLMS, WSMS, TNS, TW, SWBLR, NSPPB and NPPB. In Figure 2, group 2 is composed by characteristics related to the distance between leaves and primary branches in the main stem (ADNL and ADBPB), while group 2 in Figure 3 includes vegetative and reproductive traits related to the main stem: average distance between leaves (ADNL), number of seeds per pod on the main stem (NSPMS) and ratio between length and width of pod on the main stem (RBLWPMS). Only one characteristic (PSMS) defines group 3 in both figures. In Figure 2, group 4 includes characteristics related with total number of seeds (TNS), proportion of leaves with branches (PBL) and percentage of seed on the primary branches (PSPB). In Figure 3, group 4 includes vegetative and reproductive traits related to the main stem (NSMS, NPMS, NPB and NLMS). Group 5 (Figure 2), includes traits related to the reproductive capacity on the main stem (TSWMS, TSWPB, PWMS and PLMS). The same group on Figure 3 is related with pod and seed size (thousand seeds weight, pod and seed size parameters) and traits that include primary branches (ALPB, SLPB and WSPPB). Group 6 is mostly defined by vegetative characteristics (TBL, ALPB, HUFF and DUF).
Cluster I represent three JKI accessions in both figures, which is discriminated by morphological groups 5 and 6 (four reproductive and five vegetative characteristics) in Figure 2. In Figure 3 Cluster I is discriminated by morphological groups 1 (seven reproductive characteristics), 4 (with three reproductive and two vegetative characteristics) and 5 (five reproductive and two vegetative characteristics). Cluster II ( Figure 2) is composed by 6 accessions defined by groups 1 and 4 (five reproductive and one vegetative characteristics) and, in Figure 3, by 19 accessions and does not exhibit a defined pattern. In both figures, Cluster III with accession LM268 only, is characterized by high values in groups 1 and 5 (SWBLR, TW, PWMS, TSWPB and TSWMS) in Figure 2 and high values in group 5 and low values in groups 2-4 in Figure 3. Cluster IV in Figure 2 does not exhibit a defined pattern. For Figure 3, Cluster IV is represented by the L. albus cultivars and is characterized by high values for all heat map in most characteristics. This cluster is defined by three groups of morphological traits: group 1 (seven reproductive variables), group 2 (three reproductive and one vegetative characteristics) and group 5 (five reproductive and two vegetative characteristics).    Principal component analysis (PCA) confirmed the cluster analysis results (Figure 2 and 3). For instance: cluster I is localized oppositely for many vectors of groups 5 and 6 (TSWMS, TSWPB, PWMS PLMS, TBL, ALPB, NLMS, HUFF, and DUF) for Figure 4 and a similar scenario can be observed in Figure 5 were the vectors defined by groups 2 and 4 (HUFF, PWMS, ALPB, NPB, NLMS, SLPB, Principal component analysis (PCA) confirmed the cluster analysis results (Figures 2 and 3). For instance: cluster I is localized oppositely for many vectors of groups 5 and 6 (TSWMS, TSWPB, PWMS PLMS, TBL, ALPB, NLMS, HUFF, and DUF) for Figure 4 and a similar scenario can be observed in Figure 5

Diversity Assessed by Molecular Markers
The six selected ISSR primers used for analysis of 23 accessions resulted in the production of 37 reproducible bands (Table 5 and Figure 6). Of those, 11 (29.7%) bands were polymorphic and the remaining 26 (70.3%) were monomorphic. The total number of bands per primer ranged between four (GT 8 YC) and eight (HVH(TG) 7 ), while the percentage of polymorphic bands per primer ranged from 0 to 50%. The average for each primer was 6.2 bands. Polymorphism information content (PIC), which is used in genetics as a measure of polymorphism for a marker locus, ranged from 0.23 (HVH(TG) 7 ) to 0.72 (AG 8 YT). Effective multiplex ratio (EMR) had its minimum value with AG 8 YC (0) and maximum in GT 8 YC (2.25). The resolving power (RP) parameter used to detect the differences between a large number of genotypes ranged from 5.48 (AG 8 YT) to 13.58 (HVH(TG) 7 ). The minimum and maximum values for marker index were registered for AG 8 YC (0) and GT 8 YC (0.54) primers, respectively.

Diversity Assessed by Molecular Markers
The six selected ISSR primers used for analysis of 23 accessions resulted in the production of 37 reproducible bands (Table 5 and Figure 6). Of those, 11 (29.7%) bands were polymorphic and the remaining 26 (70.3%) were monomorphic. The total number of bands per primer ranged between four (GT8YC) and eight (HVH(TG)7), while the percentage of polymorphic bands per primer ranged    5 and MR were included in the screening test but were rejected during selection. Eight of these primers were previous tested [53].
Agronomy 2019, 10, 21 17 of 28 between a large number of genotypes ranged from 5.48 (AG8YT) to 13.58 (HVH(TG)7). The minimum and maximum values for marker index were registered for AG8YC (0) and GT8YC (0.54) primers, respectively.   The similarity matrix was used to construct a dendrogram using the UPGMA method ( Figure  7). The cophenetic correlation was 0.9058603, revealing little loss of information with transformation of similarity matrix to dendrogram. The dendrogram reveals five distinct groups. Cluster I is composed by 10 accessions of white, blue and pink flower colors and green and purple stem. Cluster II, containing 9 accessions, can be distinguished from the first group by the absence of purple stem and flower genotypes. Unlike cluster I and II, clusters III, IV and V are composed only by accessions Figure 6. Example of ISSR amplification profiles for 23 Lupinus mutabilis accessions using the primer GA 8 YT separated on a 2% agarose gel. M-NZYDNA Ladder VII marker. Numbers adjacent to accession names refer to coding used in Figure 5.
The similarity matrix was used to construct a dendrogram using the UPGMA method (Figure 7). The cophenetic correlation was 0.9058603, revealing little loss of information with transformation of similarity matrix to dendrogram. The dendrogram reveals five distinct groups. Cluster I is composed by 10 accessions of white, blue and pink flower colors and green and purple stem. Cluster II, containing 9 accessions, can be distinguished from the first group by the absence of purple stem and flower genotypes. Unlike cluster I and II, clusters III, IV and V are composed only by accessions that exhibit green stems and blue flowers. Cluster IV is represented by one accession and cluster III and V by two accessions each.

Diversity Assessed by Genomic Traits
Lupinus albus 'Misak' was validated as a DNA standard by comparison to Solanum lycopersicum 'Stupické' (Figure 8a,b) and Raphanus sativus 'Saxa' (data not shown) and estimated at 2C = 1.35 ± 0.0076 pg (1377.6 Mbp), with an average coefficient of variation of 3.47%. that exhibit green stems and blue flowers. Cluster IV is represented by one accession and cluster I and V by two accessions each.
The L. mutabilis collection was thus analysed by comparison to L. albus 'Misak', as exemplifie in Figure 8 (panels c,d). The average L. mutabilis genome size was estimated at 2C = 2.05 pg (2001 Mbp) with a 9.2% coefficient of variation, ranging from 1897.3 Mbp for accession SBP to 2083.2 Mb for accession LM34 ( Table 6). The results from a Kruskal-Wallis test performed for genome size reve significant difference between accessions (χ 2 = 94.845, Df = 23, p value = 0.000). No single accessio showed to be statistically different from all the others, rather a continuum of accessions is depicte by the homogeneous groups produced (Table 6). The L. mutabilis collection was thus analysed by comparison to L. albus 'Misak', as exemplified in Figure 8 (panels c,d). The average L. mutabilis genome size was estimated at 2C = 2.05 pg (2001.2 Mbp) with a 9.2% coefficient of variation, ranging from 1897.3 Mbp for accession SBP to 2083.2 Mbp for accession LM34 ( Table 6). The results from a Kruskal-Wallis test performed for genome size reveal significant difference between accessions (χ 2 = 94.845, Df = 23, p value = 0.000). No single accession showed to be statistically different from all the others, rather a continuum of accessions is depicted by the homogeneous groups produced (Table 6).
Genome size is an important criterion to study evolution at the intra-specific level, helping to understand conflicting pattern between morphological traits. In this study we evaluated the associations between genome size and morphological traits using Spearman correlation analysis for all 23 accessions for the two experiments. However, no single morphological trait presented strong correlation with genome size (Figure 9).

Discussion
To assess the diversity in a tarwi germplasm collection under Mediterranean conditions,

Discussion
To assess the diversity in a tarwi germplasm collection under Mediterranean conditions, phenotypic, genetic and genomic analyses were combined, studying morphological traits, ISSR markers and genome size. In general, the morphological traits used to evaluate accessions tested in the trials showed acceptable adaptability to the Mediterranean environment assuming that productivity projected is above 1.5 t/ha, achieved under rain-fed conditions unevenly distributed during the trial periods. Similar yields were previously reported in France and Spain [54]. However, the yields obtained suggest continuing breeding to achieve higher yields. The results showed also significant differences among accessions. While additional years of field trial results would certainly improve the robustness of conclusions, the results obtained are based on traits that presented coherent values between trials.
The knowledge of the correlations between different characteristics is fundamental because it allows the accomplishment of the indirect selection of the complex characteristics that are inherited quantitatively and influenced by genetic effects [43]. In this work we report positive and significant correlation between many variables. Characteristics such as total seed weight (TW) and total number of seeds (TNS) are very important and are directly related to characteristics of reproductive development. Accession LM268 presented higher values for TW and LM34, P20993 and LM18 for TNS. Therefore, these two features can work as criteria of selection in our collection for the breeding programme or to choose the most adapted. Talhinhas et al. [43] verified positive correlation between total weight with plant height, pod width, number of primary branches, proportion of seeds on the primary branches, total number of seeds and number of pods per primary branch in L. angustifolius. Georgieva and Kosev [55] also found positive correlation between thousand seed weight and plant height. Clements et al. [56] reported positive correlation between weight and plant height of 1000 seeds in L. pilosus. Heritability is a parameter widely used by breeders to genotype selection based on phenotypic expression [57]. Morphological traits exhibiting high values of heritability are chosen for the selection based on this parameter [43,58]. High heritability values enable the identification of important features to be selected for genetic breeding. Concerning our study, tarwi accessions can be selected based on the following traits: TW, SLMS, NPMS, NLMS, TNS, HUFF, PSMS, TSWMS, TSWPB, WSPPB, PWMS, SWBLR, TBL, ALPB, ADNL, SLPB, WSMS, WSPB, NPB, and NPPB. Similar results were noticed by Talhinhas [53] for SLMS, NLMS, HUFF, PSMS, TNS, and TW in L. albus, L. angustifolius and L. luteus. Our results also corroborate those by Georgieva and Kosev [55], who found high values of heritability for pod length and total number of seeds in L. albus and L. luteus.
In the present study several accessions stood out due to their superior performance in various traits. Along with accession LM268, accessions LM18, LM27, P20993, Potosi-ISA, PRT79 and I82 were the most efficient in converting vegetative growth to seed production, although lagging behind the performance of L. albus. Accession LM268 was the only tarwi accession to produce more yield on the primary branches than on the main stem, following the pattern of L. albus cultivars.
An important result worth highlighting is that most of the L. mutabilis accessions studied concentrate their production on the main stem. This characteristic is very important because it allows adaptability of tarwi to poor growth (soil and/or climatic) conditions. This characteristic may also prove useful to avoid the indeterminate growth habit of tarwi, particularly problematic in areas without summer drought. Breeding programmes should be directed for improving levels of production on the main stem and primary branches for good soil/climate conditions but with summer drought and focus on more determinate growth plants (those concentrating production on the main stem) both for marginal areas and for areas without summer drought. To the latter, accessions such as JKI-L295 and JKI-L210 stand out, as they produced over 80% of their yield on the main stem while attaining relatively high yields (ca. 10 g per plant).
The use of molecular markers in genetic diversity studies at the intra and inter specific levels proved useful in a wide range of species [59,60]. In this study we assessed the efficiency of ISSR markers for the characterisation of genetic diversity of L. mutabilis accessions. This technique is important because it allows to make a broad screening of a collection. SSR markers are not optimised yet for L. mutabilis and the transfer of such markers from other Fabaceae to tarwi did not prove successful [61], leaving ISSRs as a valid tool for preliminary screening of germplasm collections. All six primers used in this investigation revealed a polymorphism of 30.55% for all 23 accessions. Bussell et al. [62] establish 20% as minimum of monomorphic band percentage for genetic diversity study and our study reveals 69.44% monomorphic bands. Similar results were reported by Chirinos-Arias et al. [37] assessing genetic variability among 30 accessions of L. mutabilis using eight ISSR markers, finding a total polymorphism of 58.82%. The high level of polymorphism obtained in our study is in accordance with those authors. The parameters PIC, EMR, MI and RP were used to evaluate the efficiency of ISSR primers. However, to the best of our knowledge, there are no studies on L. mutabilis assessing the effective multiplex ratio, polymorphic information content, marker index and resolving power. Results show high probability in detecting polymorphism PIC (0.72), for the primer HVH(TG 7 ). The AG 8 YG primer stood out as presenting a high RP value (13.58), being more qualified to distinguish accessions. The highest value of EMR (2.25) was obtained with primer AG 8 YC, revealing this to be most efficient. The primer AG 8 YC proved to be the most useful because it presented the highest value of MI (0.54). Several studies have been undertaken based on these techniques for selecting efficient ISSR primers in different species [63][64][65][66][67].
The 23 tarwi accessions were divided in five main genetic groups using cluster analysis by the UPGMA method ( Figure 5). However, morphological characteristics such as stem and flower colour did not exhibit regular relationships in different clusters. The existence of several distinct groups that aggregate different stem and flower colours probably reflects few differences on the genetic constitution of the accessions. On the other hand, the distinct groups can reflect into distinct morphological characteristics and variations. Talhinhas [53] suggested that low intra-specific diversity in tarwi can be related to the fact that all the accessions originated from a limited number of landraces, reflecting the recent domestication genetic bottleneck effect that is estimated to have occurred no later than 2600 years before the present time in L. mutabilis [16]. Similar result was found by Chen et al. [68] in the research done on the 105 genotypes on Vigna unguiculata. In this work we verify that the genetic variability is not correlated with phenotypic variability, indicating the need for incorporation of more molecular markers. Similar results have been reported in other species. Previous studies performed by Galek et al. [14] also did not find a relation between genetic and morphological variability in accessions of L. mutabilis. In a study aiming to evaluate genetic diversity of Nelumbo using analyses of Randomly Amplified Polymorphic DNA (RAPD) and ISSR markers, Li et al. [69] found low correlation between molecular and morphological data. Talhinhas et al. [44] assessing genetic diversity in Lupinus luteus using ISSR and Amplified Fragment Length Polymorphism (AFLP) markers did not find any correlation between morphological and molecular data.
In this work we report the existence of significant differences in the intraspecific genome size (GS) variability in 23 accessions of Lupinus mutabilis. Our results reveal that the GS ranged from 1.94 pg/2C to 2.13 pg/2C. Naganowska et al. [70], also employing flow cytometry to analyse propidium iodide-stained nuclei, evaluated the nuclear DNA content variation in the genus Lupinus and found 1.90 pg/2C for Lupinus mutabilis, although a single accession was used in that study. To the best of our knowledge, our study is the first L. mutabilis genome size intra-specific analysis, depicting an overall average size of 2.05 pg (2001.2 Mbp). Several studies have reported intraspecific differences in genome size in various species such as Glycine max, Linum austriacum and Zea mays [33,71,72]. The intraspecific variation in genome size can result from repetitive/non-coding regions, hence increasing or decreasing in satellite DNA transposable elements and ribosomal genes [73]. There are studies pointing that transposable elements are largely responsible for notable differences in genome sizes. For instance, in maize, transposable elements are responsible for 85% of differences [74]. According to Petrov [75] these elements have potential of multiplicity of 0.1-1 Mbp in a single generation. The satellite DNA can also contribute greatly to genome size differences [76]. Meanwhile, Garrido-Ramos [74] refer that genomic content variation in plants which are affected by satellite DNA can range from 0.1% to 36%. Small variation of 3.5% in nuclear DNA have been associated with ribosomal genes [77]. The maximum variation of nuclear DNA content obtained in the present research was 9.2%, a value much higher than the 2% maximum genome size variability reported for soybean [71] but smaller than the 36% variation reported for maize [72]. In light of this discussion, one may discard the possibility that differences in L. mutabilis genome size are caused by the transposable elements. Only a detailed study could unravel whether this variation is due to repeated sequence differences in satellite DNA or ribosomal genes.
Data on 37 morphological traits and genome size measurement were plotted and no correlation was observed. This is not a surprise, as similar results were also reported from other studies. For instance, Oney and Tabur [31] did not find correlation between genome size and morphological traits on the Brachypodium distachyon collected in different locations in Turkey. Realini et al. [72] observed weak association between genome size and morphological traits in maize. Recently Basak et al. [78] assessing the variation of morphological traits with the genome size in turnip found no correlation. This lack of association between morphological traits and genome size suggests that other factors are determinant on the control of such characteristics, reinforcing the view that genome size variations are mainly related to non-coding regions [79].

Conclusions
The agronomic performance of L. mutabilis in Portuguese conditions was good, assuming that the assay was conducted under rain-fed conditions. Our results highlight the accession LM268 with larger seeds and a total thousand seeds weight similar to L. albus, while also achieving the highest yield and being the only tarwi accession producing more on the primary branches than on the main stem. While high yields in lupins depend on the capacity of the plants to produce large amounts of pods and seeds on lateral branches, the indeterminate growth habit of tarwi can be undesirable, either in areas without summer drought or, on the contrary, in areas with limited growing periods where further vegetative growth may impair pod filling. To this end, JKI-L295 accession present high yield concentrated on the main stem, suggesting a semi-determinate development pattern. In either case, this accession is a key point for continued breeding. In fact, the present study has shown that tarwi is still behind white lupin in terms of its adaptability to Mediterranean conditions, namely concerning yield. The genetic diversity revealed in this study, however, prompts further breeding opportunities. Molecular marker and genome size analyses have revealed important levels of genetic/genomic diversity, which could not be related to phenotypic/morphologic diversity. This illustrates a scenario of recent domestication in the absence of a gene flow to wild relatives suggesting, however, that further exploitation of genetic diversity in this tarwi collection is possible and may provide additional sources of useful agronomic traits.