A total of 166 wild macadamia accessions representing all four species [M. integrifolia (n = 49), M. tetraphylla (n = 56), M. ternifolia (n = 23), and M. jansenii (n = 23)] and one related rainforest species from the Proteaceae, Lasjia whelanii, were selected for sequencing. Within the macadamia populations, 161 wild accessions, which were collected previously from multiple locations across the natural distribution of the four species, were grown and maintained at Nambour arboretum and ex situ germplasm centers at Nambour and Tiaro in Queensland and Alstonville in NSW (Hardner et al., 2004) and five from a private collection at Limpinwood, NSW ( Supplementary Table 1 ) in Australia. Fully expanded young macadamia leaves, of accession within these ex situ collection sites, were collected in perforated labeled cellophane bags and immediately placed under dry ice until stored in a −80°C freezer at The University of Queensland, Brisbane, Australia.

Frozen leaves were coarse pulverized under liquid nitrogen using a mortar and pestle and further fine pulverized under cryogenic conditions using a Qiagen tissue lyser (MM400, Retsch, Germany). A modified version of the cetyltrimethylammonium bromide (CTAB) DNA extraction protocol described by Furtado (2014) was used to extract genomic DNA. The quality and quantity of the DNA samples were evaluated using a Nanodrop spectrophotometer (Nanodrop Technologies, Wilmington, DE, USA) by recording the absorbance ratios at 260/280 and 260/230 followed by running a 0.7% agarose gel with SYBR safe staining (Thermo Fisher Scientific). Whole-genome short-read sequencing was undertaken by BGI Hong Kong. A PCR-free library was generated and sequencing at 150-bp paired-end reads was undertaken on the DNBSEQ-G400 sequencing platform from MGI (MGI Tech Co., Ltd, Shenzhen, China) at an expected data yield/sample of at least 25× genome size (800 Mb genome)

All sequence data were analyzed in CLC Genomics Workbench 23.0.05 (CLC-GWB, CLC-Bio, QIAGEN, Denmark, http://www.clcbio.com) using the short-read pipeline. Quality control (QC) was performed for all short-read data. Reads were trimmed using a quality score limit of 0.01 with default parameters (more than 98% of the resulted trimmed reads had a Phred score >25). A subset of quality trimmed short reads (2–13 GB) was used for chloroplast genome assembly. All chloroplast genomes were assembled using the GetOrganelle toolkit (Jin et al., 2020) exploiting SPAdes v.3.15.3, Bowtie2 v.2.4.5, and Blast v.2.11.0 as dependencies. The correct configuration of the chloroplast genome was selected with respect to the M. integrifolia (sequence: NC_025288.1) (Nock et al., 2014) available at the National Center for Biotechnology Information (NCBI) (http://www.ncbi.nlm.nih.gov/) using clone manager software (Sci Ed, USA).

The chloroplast genomes were annotated using the GeSeq online tool (https://chlorobox.mpimp-golm.mpg.de/geseq.html) with M. integrifolia (sequence: NC_025288.1) as the genome (Tillich et al., 2017). Chloroplast genomes were annotated with the following settings: Annotation options: Annotate plastid inverted repeat (IR), Annotated plastid trans spliced rps12, Annotation support: Support annotation by Chloe, Annotation revision: Keep best annotation only, BLAT search–protein search identity: 25, rRNA, tRNA, and DNA search identity: 85, HMMER profile: CDS+rRNA, ARAGORN v1.2.38- Genetic code: Bacterial/plant chloroplast, Max intron length: 3,000, tRNAscan-se v2.0.7- sequence source: organellar tRNAs, MPI-MP reference set: chloroplast land plants (CDS + rRNA) and Chloe v0.1.0- annotate- CDS + tRNA + rRNA. All genomes were imported to Geneious 2023.2.1 software (Biomatters Ltd, USA) to determine the number of genes, CDS, transfer RNAs (tRNAs), and ribosomal RNAs (rRNAs) in each sample.

Previously published annotated sequences of the nuclear genome of M. integrifolia (Nock et al., 2020) were selected as reference sequences to generate accession-specific consensus CDS of nuclear genes. CDS of M. integrifolia (GCF 013358625.1) were downloaded and imported to CLC-GWB to generate a local Blast database. The CDS of 106 Arabidopsis thaliana single-copy genes identified by Li et al. (2017) was subjected to tblastn against the M. integrifolia CDS database. We selected 56 tblastn hits with a single M. integrifolia CDS matching an A. thaliana CDS as these hits represented single-copy genes in the M. integrifolia genome. From these selected single hits, corresponding M. integrifolia CDS sequences were extracted and used as a reference to extract consensus sequences from each of the macadamia accessions and from L. whelanii. BLAST analysis using the 56 extracted M. integrifolia CDS and the CDS sequences of L. whelanii as a database resulted in the selection of 53 L. whelanii CDS as single hits that represented single-copy genes in L. whelanii. Corresponding (same as selected) 53 CDS sequences from M. integrifolia and from L. whelanii were selected for further analysis. The 53 CDS from M. integrifolia ( Supplementary Table 2 ) were used as reference sequences in the read mapping approach to generate corresponding CDS consensus sequences for each of the macadamia accessions. Essentially, short reads trimmed data of each macadamia accession were mapped separately to each of the 53 selected M. integrifolia CDS sequences to extract consensus sequences. Finally, consensus CDS sequences, were extracted for each macadamia accession, and were concatenated in the same sequential order to obtain the final nuclear gene CDS sequence.

Phylogenetic evaluation of macadamia was conducted by utilizing complete chloroplast genome sequences and single-copy concatenated nuclear gene CDS sequences. For phylogenetic analysis, we selected 138 wild macadamia accessions representing all four species. To maintain precision and clarity within our analysis, we exclude accessions from planted wild germplasm, accessions known to be natural hybrids, admixtures, accessions of unknown origin, and unidentified macadamia accessions ( Supplementary Table 3 ). The phylogenetic trees were also generated for four species: M. integrifolia (n = 44), M. tetraphylla (n = 49), M. ternifolia (n = 22), and M. jansenii (n = 23) separately based on chloroplast genomic data and single-copy nuclear gene sequences. L. whelanii was used as an outgroup. All selected sequences were aligned along with the outgroup using MAFFT alignment with default parameters in Geneious 2023.2.1 software (Biomatters Ltd, USA).

Geographical maps of origin were created for all four macadamia species based on the chloroplast phylogenetic clade separation. Geographic ranges were mapped using Esri National Geographic in the QGIS Geographic Information System (Version 3.32.1-Lima) (http://qgis.osgeo.org). Geographical coordinates of each accession are listed in Supplementary Table 3 .

Following paired-end sequencing (150 bp), a total of 140,789,058–260,592,896 reads were obtained for 166 Macadamia accessions ( Supplementary Table 5 ) with sequence depth of 18.23× and 50×. The sequence depth of trimmed paired-end reads at a quality score limit of 0.01 ranged between 16.47× and 43.98×. Two assembled sequences were obtained from GetOrganelle analysis indicating the presence of two structural haplotypes of the chloroplast genome that occurs in plants related to the orientation of the single-copy region. The correct configuration of the chloroplast genome was selected with respect to M. integrifolia (Reference sequence: NC_025288.1). Complete chloroplast genome sizes analyzed in this study are shown in Supplementary Table 6 . Complete circular chloroplast genomes were obtained for all genotypes, ranging in size from 159,195 to 159,734 bp. The smallest chloroplast genome was identified for three M. tetraphylla accessions (Mac_297, Mac_338, and Mac_345) and two wild macadamia trees of uncertain species (Mac_047 and Mac_329) while the largest was observed for two M. integrifolia accessions (Mac_029 and Mac_262). Chloroplast genome sizes of M. integrifolia ranged from 159,458 to 159,734 bp, those of M. ternifolia ranged from 159,463 to 159,508 bp, those of M. tetraphylla ranged from 159,195 to 159,598 bp, and M. jansenii were 159,524 bp in length except for MacP_16 (159,526 bp). The result of this study also revealed that out of 23 M. jansenii accessions, 22 accessions had identical chloroplast genomes.

All macadamia chloroplast genomes showed a quadripartite structure of angiosperm, including a large single copy (LSC), a small single copy (SSC), and two identical inverted repeats (IRa and IRb) ( Figure 1 ). Gene annotation showed 116 full-length genes, 81 CDS, 4 rRNAs, and 31 tRNAs. Among these 116 genes, 60 genes were involved in protein synthesis and DNA replication (genes responsible for rRNAs, tRNAs, large subunit of ribosome, small subunit of ribosome, and DNA-dependent RNA polymerase), 46 were involved in photosynthesis (genes responsible for subunits of photosystem I, subunits of photosystem II, subunits of ATP synthase, subunits of NADH dehydrogenase, large subunit of rubisco, and subunits of cytochrome complex), 6 were involved in other different functions (genes responsible for inner membrane protein, cytochrome synthesis gene, acetyl-CoA-carboxylase, maturase, ATP-dependent protease, and translational initiation factor), and 4 were involved in unknown function genes ( Supplementary Table 7 ).

The L. whelanii chloroplast genome possessed the standard quadripartite structure, containing two inverted repeats (18,824 bp), the LSC region (87,911 bp), and the SSC region (26,448 bp). ( Figure 1 ). Genome size was recorded as 159,631 bp. The plastome of L. whelanii contained no significant difference in relation to genes, protein coding genes, rRNA, and tRNA. Overall, the GC content of the chloroplast genome was recorded as 38%.

The multiple chloroplast genome alignment of 44 M. integrifolia accessions together with the outgroup L. whelanii was 161,281 bp in length with 99.6% identical sites. The tree topologies of both were similar ( Figure 2A ), and most nodes were supported by high bootstrap support (BS) (>95%) and Bayesian posterior probabilities (PP) (>0.95). However, some internal nodes tended to have low BS, indicating incomplete lineage sorting. The phylogenetic tree construction revealed that Mac_232 (corresponding to population site 90) clusters separately from the rest of the 43 accessions. The remaining accessions were differentiated into two main clades and further differentiated into sub-clades. Clade II contained accessions from the northern distribution of M. integrifolia: Mac_231, Mac_262, Mac_029, Mac_265, and Mac_033 from the Gundiah/Mount Bauple region (corresponding to population sites 1, 2, 2, 3, and 3, respectively) ( Figure 2B ) and Mac_052, Mac_091, Mac_340, Mac_248, and Mac_266 from the Gympie region (corresponding to population sites 9, 55, 56, 57, and 57, respectively). Clade III included a total of nine accessions, of which seven were from the Caboolture region: Mac_250, Mac_312, Mac_246, Mac_045, Mac_080, Mac_026, and Mac_143 (corresponding to population sites 20, 21, 21, 71, 76, 76, and 77, respectively) and one from the Nambour region: Mac_044 (corresponding to population site 101). Interestingly, Mac_059 from population site 57 did not cluster with two other accessions (Mac_248 and Mac_266) from population site 57 (Gympie region) in Clade I. Clade IV contained 24 accessions that were collected from the region south of Brisbane except for Mac_251, Mac_089, Mac_139, and Mac_189 (corresponding to population sites 3, 90, 90, and 90, respectively). This observation verified that the accessions having the same geographical origin tend to form distinct clusters among themselves. However, it is noteworthy that certain accessions from these localities exhibit a tendency to associate or cluster with accessions sourced from different geographical areas. The result of this study also revealed that Mac_251, Mac_089, Mac_139, and Mac_189 (corresponding to population sites 3, 90, 90, and 90, respectively) are likely to be planted trees. Mac_312, Mac_246 and Mac_306, Mac_228 are biological replicates and confirmed by clustering together. Additionally, the short branch length of the tree suggested that Mac_232 had less divergence while Mac_026, Mac_044, Mac_080, and Mac_143 from Clade III were more diverged accessions.

To study the phylogenetic relationship of M. tetraphylla, phylogenetic trees were constructed using 49 complete chloroplast genome sequences. The multiple chloroplast genome alignment of M. tetraphylla accessions together with the outgroup L. whelanii was 161,554 bp in length with 96% identical sites. ML and BI trees exhibited similar phylogenetic topologies ( Figure 3A ). The resulting phylogenetic tree showed strong statistical support for most internal and external nodes but barring poor BS value for some internal nodes. However, in the Bayesian chloroplast tree, the highest PP value of 1 was observed for all the nodes. The chloroplast tree displayed two major clades. The first major clade (Clade I) consisted of germplasm collected from the southern part: the Lismore region (Mac_060, Mac_108, Mac_134, Mac_268, Mac_244, Mac_031, Mac_297, Mac_184, Mac_083, and Mac_095) and the Ballina region (Mac_247, Mac_325, MacP_14, and Mac_115) except for MacP_15 (corresponding to population site 31) and Mac_097 (corresponding to population site 38) ( Figure 3B ). The second major clade was further divided into sub-clades. Clade II contained accessions from the Murwillumbah region: Mac_341, Mac_314, Mac_227, Mac_270, Mac_098, Mac_064, Mac_291, and Mac_259 (corresponding to population sites 37, 37, 81, 160, 160, 160, 160, and 160, respectively) and the Beenleigh region: Mac_236 (corresponding to population site 100). Two accessions from Clade III from population site 37 clustered separately from the rest of the accessions from the same geographical location. However, as in the M. integrifolia chloroplast phylogenetic tree, the majority of M. tetraphylla tended to form distinct clusters among themselves based on geographical areas. Results also indicated that Mac_031, Mac_244, and Mac_268 from population site 96 (Clade I) and Mac_264 and Mac_238 from population site 84 (Clade VI) were highly diverged accessions.

The phylogenetic relationships within M. ternifolia accessions were inferred by 22 assembled complete chloroplast genomes. A multiple chloroplast alignment conducted using an outgroup was 160,807 bp with 97.2% identical sites. Phylogenetic trees built with the whole chloroplast genome using both methods had the same topology ( Figure 4A ). The results showed two major clades having MacP_11, Mac_309, and MacP_12 from Nambour (corresponding to population site 88) ( Figure 4B ) in one major clade and the remaining accessions in the second major clade. There was a clear relationship between the phylogenetic structure and geographic origin of the wild accessions of M. ternifolia. The resulting topology suggested accessions from Clade IV: Mac_299, Mac_071, Mac_317, Mac_332, and Mac_334 (corresponding to population sites 20, 51, 51, 51, and 51, respectively) were highly diverged accessions. In the M. jansenii population, no variants were observed, indicating that all accessions shared the same chloroplast haplotypes except for MacP_16 with a 2-bp difference.

The chloroplast phylogeny tree generated by taking 138 complete chloroplast genomes was supported with a PP of 1.0. Two major clades were identified ( Figure 5 ). The first major clade contained 16 M. tetraphylla accessions from the Lismore region (Mac_060, Mac_108, Mac_134, Mac_268, Mac_244, Mac_031, Mac_297, Mac_184, Mac_083, and Mac_095), the Ballina region (Mac_247, Mac_325, MacP_14, and Mac_115), the Beenleigh region (MacP_15), and the Murwillumbah region (Mac_097). Interestingly, all these accessions corresponded to Clade I in M. tetraphylla chloroplast phylogenetic tree ( Figure 3A ). The second major clade was further differentiated into two sub-clades. All the M. jansenii accessions were clustered in one clade. The second sub-clade was further divided into two clades. The small sub-clade contained 10 accessions from the northern distribution of M. integrifolia (corresponding to Clade II in the M. integrifolia Cp phylogenetic tree) and 1 accession from the Nambour region (corresponding to Clade I in the M. integrifolia Cp phylogenetic tree). The larger sub-clade contained all the remaining M. tetraphylla, M. integrifolia, and M. ternifolia. This result shows that accessions that were collected from same locality cluster together. We assumed that chloroplast capture could be the reason for the presence of different species in the same clade when a species coexists in the same geographic area with other species.