Lima bean (*Phaseolus lunatus* L.) is one of five domesticated *Phaseolus* species, and alongside common bean, holds significant agronomic and economic importance. Its wide ecological adaptability across diverse climates from Mexico to Argentina makes it a promising crop for enhancing food security, especially under projected climate change scenarios in Latin America and beyond. Lima bean seeds are rich in protein and carbohydrates, containing essential amino acids. Wild and domesticated Lima beans show significant diversity, structured into three major gene pools (two Mesoamerican – MI and MII – and one Andean – AI). Evidence suggests at least two independent domestication events, one in Mesoamerica and another in the Andes. Both domestication events resulted in convergent evolution, leading to similar traits such as larger pods and seeds, reduced pod dehiscence, and loss of seed dormancy. Previous research relied heavily on common bean reference genomes, which can introduce bias. Therefore, a dedicated Lima bean reference genome is crucial to avoid reference bias and unlock deeper insights into its genetic diversity and domestication.

Previous research has established the existence of three major gene pools in wild Lima beans: two Mesoamerican (MI and MII) and one Andean (AI). These gene pools have mostly non-overlapping distributions. Studies using genomic data have supported the hypothesis of at least two independent domestication events, one in Mesoamerica and one in the Andes. The Andean domestication event involved gene pool AI, resulting in varieties with large, flat seeds (Big Lima cultivars), while the Mesoamerican domestication event arose from gene pool MI, producing varieties with smaller, rounded or oval seeds (Potato and Sieva cultivars). The convergent evolution of similar traits in both Andean and Mesoamerican landraces highlights the species’ adaptability. The high levels of homozygosity found in both Lima and common beans due to predominantly self-pollinating nature has made the creation of high quality reference genomes more feasible, although the lack of previous genomic resources meant researchers had previously relied upon common bean genomes. This created a research bias.

This study employed a multi-faceted approach combining PacBio and Illumina sequencing technologies to generate a high-quality, chromosome-level genome assembly for a domesticated Lima bean accession (G27455) from the Mesoamerican gene pool MI. Four sequencing protocols were utilized: paired-end whole-genome sequencing (WGS), 10x Genomics linked-read sequencing, genotyping-by-sequencing (GBS), and RNA sequencing. A linkage map was constructed using GBS data from an F5 generation of a recombinant inbred line (RIL) population (UC 92 × UC Haskell). This linkage map enabled scaffolding and ordering of contigs into 11 pseudomolecules representing the Lima bean chromosomes. Gene models were predicted using a combination of ab initio prediction, RNA-Seq data from leaf, flower, and pod tissues at different developmental stages (two accessions), and homology-based approaches. Repetitive elements were masked, and functional annotation was performed using various tools (RepeatMasker, Maker, Trinotate). Comparative genomics was conducted with the common bean genome to identify orthologs, paralogs, and structural variations. QTL mapping was performed on the RIL population to identify genomic regions associated with traits such as determinacy, flowering time, hundred-seed weight, and cyanogenesis. Population structure analysis was conducted using GBS data from nearly 500 wild and domesticated accessions. Different clustering methods (Neighbor Joining, Discriminant Analysis of Principal Components, STRUCTURE, fineSTRUCTURE) were employed to infer population structure and relationships among accessions. RNA-seq was performed to compare gene expression in pod development between a wild and domesticated accession at two stages (T1: initiation of pod elongation, T2: before seed filling) to investigate the genetic basis of pod dehiscence. Functional enrichment analysis was done to identify pathways associated with differentially expressed genes.

A high-quality chromosome-level genome assembly (542 Mbp) of Lima bean was generated. 28,326 gene models were predicted and annotated, with 1917 genes identified as potentially involved in disease resistance, many of which showed high diversity. Comparative genomics with common bean revealed high synteny but also five major intrachromosomal rearrangements. Population genomic analyses identified six distinct clusters within wild Lima bean, mainly with non-overlapping distributions. Mesoamerican landraces were further subdivided into three sub-clusters. RNA-seq revealed 4275 differentially expressed genes (DEGs) related to pod dehiscence and seed development. QTL mapping identified nine QTLs associated with agronomic traits, including a major QTL for determinacy on chromosome P101 (likely *PITFL1*), a major QTL for flowering time on P101 (*PITFL1*), four minor QTLs for hundred-seed weight and three QTLs for cyanogenesis on chromosomes P105, P108, and P110 showing epistatic interactions. Linkage disequilibrium (LD) decayed faster in wild than domesticated accessions. Introgression analysis revealed instances of chromosomal segment exchange between gene pools, potentially due to recent contact between wild and domesticated populations or between domesticated populations of different origins. Analysis of gene expression during pod development identified a significant increase in expression of the *PDH1* gene (involved in lignin deposition) at T2 in the wild accession, consistent with its role in pod dehiscence.

The findings address the research question by providing a comprehensive genomic resource for Lima bean. The high-quality genome assembly and extensive annotation are significant contributions to the field, enabling future studies on gene function and evolution. The identification of disease resistance genes and QTLs for agronomic traits provides valuable information for breeding programs. The insights gained into the population structure and domestication history of Lima bean improve our understanding of its evolutionary trajectory and the genetic basis of domestication traits. The observed convergent evolution of traits in different domestication events highlights the species' adaptability. The RNA-seq data provide functional insights into the molecular mechanisms underlying important traits such as pod dehiscence. The identification of specific genes and their regulatory patterns offers opportunities for targeted breeding strategies to improve yield and other traits. The study's comprehensive analysis illuminates the genetic architecture of this important crop, paving the way for further research and development of improved varieties.

This study provides a comprehensive set of genomic resources for Lima bean, including a high-quality genome assembly, gene annotation, QTL mapping results, and population genomic data. These resources will facilitate future research on the evolution, genetics, and breeding of Lima bean. Further research should focus on functional validation of candidate genes identified, further exploration of the genetic basis of climate change adaptation, and investigation of gene interactions to fully understand the complex architecture of important agronomic traits. The development of marker-assisted selection strategies based on the identified QTLs will contribute towards more efficient breeding for improved varieties.

The study primarily focused on a single domesticated accession for genome assembly, potentially limiting the representation of genetic diversity. The GBS data, while providing a large number of accessions, has inherent limitations in terms of SNP density and potential missing data, which might affect the accuracy of some analyses. The RNA-seq analysis focused on a limited number of tissues and developmental stages, and the generalization to other tissues and stages may require additional investigation. The use of one wild and one domesticated accession for the pod dehiscence study limits inferences about population-wide variation.