Complete Mitochondrial Genomes Reveal Population-Level Patterns in the Widespread Red Alga Gelidiella fanii (Gelidiales, Rhodophyta)

Although complete mitogenomic data have been widely applied in human and other animal population studies, they are extremely limited for florideophycean red algal populations. Gelidiella fanii is a recently described rhodophyta, previously misidentified as G. acerosa, a cosmopolitan agar-yielding species from tropical to subtropical waters. To decipher patterns in genetic diversity and geographic distribution for G. fanii, we obtained 10 complete mitogenomes including two outgroups, G. acerosa and G. flabella. The mitogenomes ranged in size from 25,223 to 25,281 bp and had 48 genes, which are similar in general structure, gene order and content, and presence of a group II intron. Phylogenomic analysis revealed that G. fanii was monophyletic and clearly separate from G. acerosa. The range of G. fanii was extended from Southeast Asia and northern Australia to Eritrea, Juan de Nova Island, and Kenya in the west, and to Hawai‘i and Tetiaroa Atoll to the east. Haplotype network analysis of cox1 revealed seven geographically structured groups: Southeast Asia, Kenya/Juan de Nova Island, Indonesia, northern Australia, the Philippines, Tetiaroa Atoll, and Hawai‘i. This regional structure has likely resulted from the separation and isolation of an ancient widespread population during the Pleistocene. This study demonstrates that mitogenome sequencing is a powerful genotyping tool for studies of genetic diversity, biogeography, and conservation of economically valuable marine algal species.


INTRODUCTION
In red algae, the mitochondrial genome (mitogenome) comprises circular, maternally inherited chromosomes with fast evolving genes (Hughey et al., 2014;Yang et al., 2015). First sequenced from Chondrus crispus Stackhouse (Leblanc et al., 1995), mitogenomes have improved our knowledge of evolution, genetics and the taxonomy of red algae. Hughey et al. (2014) was the first to analyze mitogenomes from archival type specimens of bangiophycidean Pyropia species, and proposed a genome-based metric for distinguishing species. Yang et al. (2015) analyzed the mitogenomes of representatives of all five subclasses in the Florideophycidae. They documented the rapid radiation of the class, but concluded that gene synteny in these phylogenetically diverse red algae was highly conserved. In a review of 16 mitogenomes among 30 datasets publicly available at that time, Salomaki and Lane (2016) reported that red algal mitogenomes were more highly conserved than previously reported. Boo et al. (2016a) analyzed 10 mitogenomes from type specimens of species in the order Gelidiales and proposed that, of 23 protein-coding genes (PCGs) present, six PCGs with low nonsynonymous (Ka)/synonymous (Ks) ratios were suitable markers for species identification. Mitogenomes have also been used for the purpose of merging morphological species, describing new species and genera, and identifying introduced species (Hughey and Boo, 2016;Suzuki et al., 2016;Song et al., 2017;Gabrielson et al., 2018;Boo and Hughey, 2019;Bustamante et al., 2019). However, red algal mitogenomes have yet to be analyzed from a population perspective.
Gelidiella fanii S.-M.Lin (Gelidiellaceae) is one of 15 species in the widespread genus Gelidiella Feldmann and Hamel, found on coral reefs and intertidal rocky shores in tropical to subtropical waters (Abbott, 1999;Costa et al., 2002;Boo et al., 2015;Huisman et al., 2018). Gelidiella fanii is morphologically similar to G. acerosa (Forsskål) Feldmann and Hamel, a cosmopolitan species yielding high-grade bacteriological and pharmaceutical agarose, as well as agar for food (Rioux and Turgeon, 2015). Gelidiella fanii is a recently described species from Taiwan that was previously misidentified as G. acerosa (Lin and Freshwater, 2008). It has been shown to be clearly distinct using mitochondrial cox1 and plastid rbcL DNA sequences, as well as by morphology. Thalli are iridescent under water and have downward curved branches with slender unilateral branchlets, numerous surface hairs on the distal end of branches and branchlets, and smaller tetrasporangia (Lin and Freshwater, 2008;Wiriyadamrikul et al., 2010). There have been no reports of sexual reproduction in G. fanii (Lin and Freshwater, 2008;Wiriyadamrikul et al., 2010;Boo et al., 2016c;Huisman et al., 2018), a characteristic common in other species of Gelidiella. Molecular analyses have confirmed the presence of G. fanii in Indonesia, Japan, the Philippines, Thailand, Vietnam, and northern Australia, in addition to Taiwan (Lin and Freshwater, 2008;Wiriyadamrikul et al., 2010;Boo et al., 2016c;Huisman et al., 2018). However, the taxonomy and distribution of G. fanii remains understudied outside of Southeast Asia and northern Australia because of its morphological similarity to G. acerosa, a species that is recorded globally and is morphologically variable depending on habitat.
To date, 40 individual gene sequences of G. fanii are publicly available in GenBank; 23 mitochondrial cox1, and 14 plastid rbcL, two psaA, and one psbA from Southeast Asia and northern Australia (Lin and Freshwater, 2008;Wiriyadamrikul et al., 2010;Boo et al., 2016b,c;Huisman et al., 2018). These data have been used for taxonomic studies and analyses of the biogeographical structure of cox1 haplotype networks in Southeast Asia. However, limited taxon sampling and individual gene datasets are insufficient for full molecular characterization of the populations of this presumably widespread species. We wished to determine how complete mitogenomes vary among geographically isolated populations of G. fanii, and which individual genes have useful resolution at the population level. The aims of this study were i) to obtain complete mitogenomes using High Throughput Sequencing (HTS) techniques, ii) to investigate the utility of intergenic spacer regions as well as PCGs for population studies, and iii) to establish haplotype networks of cox1 to understand the distribution of G. fanii in the Indo-Pacific Ocean.

Habitat and Collection of Specimens
Like G. acerosa, G. fanii occurred on subtidal coral reefs and/or intertidal rocky reefs in the Indo-Pacific regions; however, G. fanii was very rare compared with G. acerosa. Specimens were mostly typical of the species, but they sometimes displayed fewer unilateral branchlets. All specimens collected in this study were vegetative or tetrasporic; and neither spermatangial (male) nor cystocarpic (female) plants were found.
Fresh specimens were collected in Hawai'i, Juan de Nova Island in the Mozambique Channel, Tetiaroa Atoll in French Polynesia, the Philippines, Taiwan and Vietnam, and were placed in individual bags with silica gel until processed. Specimens were identified based on morphological observation as well as analyses using mitochondrial cox1 and plastid rbcL sequences. In addition, herbarium specimens identified as G. acerosa were studied on loan from University of California at Berkeley, United States. Voucher specimens are housed in the Natural History Museum, Chungnam National University, Korea (CNUK) and the herbarium of the University of French Polynesia, Tahiti (UPF). Information on all specimens used in this investigation is provided in Supplementary Table S1.

DNA Extraction and Genome Sequencing
Eight Gelidiella fanii individuals were selected to represent the species' broad geographic distribution in the Indo-Pacific Ocean and the topology of the cox1 phylogeny (below). DNA was extracted from ∼5 mg of dried tissue using NucleoSpin Plant II Kit (Macherey-Nagel, Düren, Germany) according to the manufacturer's protocol or DNeasy Blood and Tissue Kit (Qiagen, Valencia, California, USA) following Boo et al. (2016a). Two genomic DNAs extracted in the previous studies (Wiriyadamrikul et al., 2010;Huisman et al., 2018), one from the Philippines (CNU026931) and the other from northern Australia (CNU066493), were also used for genome sequencing. Two outgroup species, G. acerosa and G. flabella G.H.Boo and L.Le Gall, were included for HTS. Library preparation and HTS were performed by Genotech Co. (Daejeon, Korea) with Illumina platforms, HiSeq2500 or HiSeqX, using 100 bp or 150 bp paired-end library constructions, respectively.

Genome Assembly and Annotation
The raw reads were assembled using a combination of approaches, NOVOPlasty 3.5 (Dierckxsens et al., 2016), SPAdes assembler 3.13.0 (Nurk et al., 2013) and MEGAHIT . Assembled contigs were sorted and reassembled using Geneious Prime 2019.2.3 1 to construct consensus mitochondrial genome sequences. To confirm the accuracy of the assembly, raw reads were mapped to the draft mitochondrial genomes using Bowtie2 in Geneious Prime.

Comparison of Mitogenome Structure and Phylogenetic Analysis
The physical map of the mitogenome was prepared for visualization using OrganellarGenomeDRAW (OGDraw) (Lohse et al., 2013). Locally collinear blocks (LCBs) alignments were generated using ProgressiveMauve (Darling et al., 2004) with a seed of 21 for the mitochondrial alignments and the 'Use seed families' option selected. CREx (Bernt et al., 2007) was performed to compare the gene order and their rearrangement events (e.g., inversion, transpositions, reverse transpositions, and tandem-duplication-random-loss [TDRL]) using heuristic pairwise comparisons with the common interval measurement in the genus Gelidiella. Twenty-three PCGs were translated into amino acid sequences and aligned using MAFFT 7.450 (Katoh and Standley, 2013) under the default setting in Geneious Prime. The alignment was manually adjusted and the ambiguous sites were deleted by GBlocks v.0.19b (Castresana, 2000). The phylogenetic tree of the concatenated dataset of 23 PCGs was reconstructed using Maximum Likelihood (ML) and Bayesian inference (BI). For the ML analysis, the best-fitting partitioning schemes and models of molecular evolution were inferred by using ModelFinder (Kalyaanamoorthy et al., 2017). Based on the Bayesian information criterion (BIC), ModelFinder identified two partitions: mtVer + F + I + G4 for atp4 + atp8 + rpl16 + rps3 + rps11 + sdh3 + sdhD + tatC and mtVer The ML analyses were performed using the W-IQ-tree webserver (Trifinopoulos et al., 2016) with 1,000 ultrafast bootstrap (BS) replicates. For the BI analysis, PartitionFinder v.2.1.1 (Lanfear et al., 2016) was used to select the best-fitting partitioning schemes and models of molecular evolution using the greedy algorithm with unlinked branch lengths. The PartitionFinder identified three partitions: CPREV + I + G for nad2 + nad5 + nad3 + nad4 + sdh2 + rps12 + nad6 + rpl16 + atp8 + rps11 + rps3 + atp4, MTMAM + I + G for cox1 + atp9 + cob + cox2 + cox3 + atp6 + nad1 + nad4L, and MTMAM + I + G for tatC + sdh3 + sdhD. The BI was performed with MrBayes v.3.2.1 (Ronquist et al., 2012) using the Metropolis-coupled Markov Chain Monte Carlo (MC3) with the models selected by PartitionFinder. For each matrix, four million generations of two independent runs were performed with four chains and sampling trees every 100 generations. The burn-in period was identified graphically by tracking the likelihoods at each generation to determine when they reached a plateau. Twenty-five percent of saved trees were removed, and the remaining trees used to calculate the Bayesian posterior probabilities (BPP).

Molecular Analyses of Mitochondrial cox1 and Plastid rbcL
DNA extractions, polymerase chain reaction (PCR) amplification, and sequencing followed Boo et al. (2016b). The primers used for amplifying and sequencing were F7, F645, R753, and RrbcS start for rbcL (Freshwater and Rueness, 1994;Lin et al., 2001;Gavio and Fredericq, 2002), and COXI43F and COXI1549R for cox1 (Geraldino et al., 2006). When large fragments of the analyzed loci could not be amplified in one herbarium specimen (UC1461694, collected in 1962), we were able to amplify and sequence 237 bp of rbcL using primers F577, F993, and R753 (Freshwater and Rueness, 1994). Sequences of the forward and reverse strands were determined for all taxa, and the electropherograms were edited using MEGA7 (Kumar et al., 2016) and checked manually. Newly generated sequences were deposited in GenBank. Sequences were aligned using the MUSCLE algorithm in MEGA7 with default parameters and the alignment was manually adjusted.
Phylogenies of individual datasets were reconstructed with Maximum Likelihood (ML) analysis using the W-IQ-tree webserver. The best-fitting substitution model was determined with the model test option (auto), followed by the ML tree search, and 1,000 ultrafast bootstrap replicates. Haplotype network of mitochondrial cox1 was constructed with PopART v.1.7 using the median-joining networks (MJN). Haplotype diversity (h) and nucleotide diversity (π) were calculated for each population and at the species level using DnaSP v.5.1. Non-hierarchical and hierarchical analyses of molecular variance (AMOVA) was performed using Arlequin v.3.5 (Excoffier and Lischer, 2010) with -statistics to quantify the proportion of total genetic variance, with significance of fixation indices tested using 10,000 permutations. The hierarchical partition was set to seven groups (I-VII) based on the phylogeny and haplotype network of cox1.

Gelidiella Mitogenomes
Ten complete mitogenomes were sequenced using HTS; eight from Gelidiella fanii and one of each from G. acerosa and G. flabella ( Table 1). The sequences of Gelidiella mitogenomes showed 91.5% identity in most regions relative to G. fanii. Pairwise divergences of mitogenomes within G. fanii were in a range of 2.4 ± 0.8% (between Hawai'i and northern Australia). Interspecific pairwise divergences were 14.4 ± 0.1% between G. fanii and G. acerosa and 18.0 ± 0.08% between G. fanii with G. flabella.
The mitogenomes of G. fanii ranged in size from 25,223 bp (JN071, Juan de Nova Island) to 25,263 bp (CNU026931, the Philippines), with highly conserved gene synteny (Table 1 and Supplementary Figure S1). The GC content was in a range of 30.3 ± 0.1%. The mitogenomes contained 48 genes, consisting of 23 PCGs, 23 tRNAs, and 2 rRNA subunits, a result similar to the publicly available mitogenomes of Gelidium J.V.Lamouroux and Pterocladiella B.Santelices and Hommersand (Boo et al., 2016a), which were included in our phylogenomic analysis. All eight G. fanii mitogenomes had a group II intron between the nad5 and nad4 genes. Three types of tRNAs were found between the trnA and trnN regions ( Figure 1A). TYPE1 (trnY-trnR-trnS insertion) was found in most G. fanii and G. acerosa samples. TYPE2 (trnS-trnR-trnY inversion) was present in G. fanii from northern Australia (CNU066493). This tRNA inversion rearrangement was inferred from the putative ancestral state of G. acerosa (trnY-trnR-trnS). TYPE3 (lacking trnY-trnR-trnS) was most common in the Gelidiales and Gelidiella flabella ( Figure 1A and Supplementary Table S2).

Phylogenomics and Gene Characteristics
The concatenated dataset of 23 PCGs (5,555 amino acid positions) from 23 mitogenomes, including 13 previously published mitogenomes (nine Gelidium and four Pterocladiella), was used for phylogenetic analysis. Because the topologies of ML and BI were identical, we show the ML tree with branch supports of MLBS and BPP. Mitochondrial phylogenomics revealed that the three species of Gelidiella formed a fully supported clade. Gelidiella fanii was monophyletic (BS: 100, BPP: 1.0) and clearly segregated from G. acerosa and G. flabella (Figure 1B). Two G. fanii mitogenomes from the Philippines clustered with full support with those from Tetiaroa Atoll (BS: 100, BPP: 1.0) and those formed a moderately supported clade with Hawai'i (BS: 79, BPP: 0.97). Gelidiella fanii from Taiwan formed a clade (BS: 100, BPP: 1.0) with Vietnam, and was sister to Juan de Nova Island.
The characteristics of the 23 PCGs are provided in Table 2. The mean value of the ratio of non-synonymous (Ka) versus synonymous substitutions (Ks) for 23 PCGs was in a range of 0.0000-0.2827 (Supplementary Figure S2). The value was low in atp9 (0.0000), cox1 (0.0189), cob (0.0258), cox2 (0.0416), cox3 (0.0519), and nad4L (0.0538) compared to the other PCGs. The number of haplotypes for 23 PCGs is provided in Table 2. The haplotype networks of six genes with a low Ka/Ks ratio are shown in Supplementary Figure S3.

Phylogeny of Mitochondrial cox1 and Plastid rbcL
Both markers newly confirmed the occurrence of G. fanii in Kenya, Juan de Nova Island, and Tetiaroa Atoll; rbcL sequences confirmed its presence in Eritrea (Supplementary Figure S4). The unrooted phylogeny of G. fanii based on 47 cox1 sequences, including 21 newly generated sequences, consisted of seven groups (I-VII) (Figure 3); I from Southeast Asia (Taiwan, Thailand, the Philippines, Vietnam, Japan), II from Eastern Africa (Juan de Nova Island, Kenya), III  from Indonesia, IV from northern Australia, V from the Philippines, VI from Tetiaroa Atoll, and VII from Hawai'i. Groups III and IV formed a strongly supported clade (BS: 94); groups V and VI were clustered with moderate support (BS: 84). Of a total of 23 rbcL sequences from G. fanii, 10 were generated in this study. The short rbcL sequences (237 bp) were successfully amplified from the herbarium specimen from Eritrea (UC1461694, as G. acerosa) that we identified as G. fanii. It differed by 0 to 1 bp (0-0.4%) from all other G. fanii specimens (Supplementary Table S3), but was 6 bp (2.5%) different from G. acerosa and 16 bp (6.8%) from G. flabella. The unrooted phylogeny of rbcL was similar to that of cox1 (Supplementary Figure S4), except that the Hawaiian sequences grouped with Kenya and Juan de Nova Island.
Genetic Diversity, Haplotype Network, and Population Structure of Mitochondrial cox1 The cox1 dataset revealed high estimates of haplotype diversity (h = 0.931) and nucleotide diversity (π = 0.00951), indicating genetic heterogeneity within the species. The median-joining network consisted of 19 haplotypes in seven groups, connected by many missing haplotypes related to the geographical distances between sampling locations (Figure 4). Haplotypes from Indonesia (C13) and northern Australia (C14) were distantly related (32 missing haplotypes) while Southeast Asian haplotypes (C1-C10) were connected by 1-3 missing haplotypes. Two haplotypes from Taiwan shared sequences with Thailand (C2) and the Philippines (C6), respectively. Haplotypes from the Philippines were separated into two groups: one was restricted to Cebu (C15-C17), while the other (C6), shared with Taiwan, clustered with Southeast Asian haplotypes. The Indonesian haplotype (C13) was distantly connected to those from Southeast Asia. Haplotypes from northern Australia (C14), Tetiaroa Atoll (C18), and Hawai'i (C19) were also isolated. Non-hierarchical AMOVA showed that most of the cox1 variation within species was found among populations (85.36%; Table 4), while a smaller amount of genetic variation was found within populations (14.64%). Genetic subdivision was highly significant among populations ( ST = 0.85, P < 0.001). The hierarchical AMOVA indicated the total genetic variance was mainly explained by the variance among groups (88.16%, CT = 0.88, P < 0.001); the variances among populations within groups (6.23%, SC = 0.53, P < 0.001) and within populations (5.62%, ST = 0.94, P < 0.001) were much lower ( Table 4).

DISCUSSION
We sequenced eight complete mitogenomes of Gelidiella fanii and two outgroup species, G. acerosa and G. flabella. Gene content and organization in the newly determined Gelidiella mitogenomes are similar to published reports for Gelidium and Pterocladiella. However, the Gelidiella mitogenome is slightly larger (45-48 genes and 25,223-25,281 bp) than those (43-44 genes and 24,901-24,970 bp) of Gelidium and Pterocladiella (Yang et al., 2015;Boo et al., 2016a;Boo and Hughey, 2019). The differences are due to the addition of three tRNA genes centrally positioned in the tRNA track between atp6-tatC of all three species of Gelidiella. The three tRNAs are tyrosine, arginine, and serine, which are located between arginine and asparagine. Interestingly, an insertion of tRNA track, trnY-trnR-trnS (TYPE1) was characteristic of G. fanii and G. acerosa, supporting their close relationship compared to G. flabella and other species in the Gelidiales. The Australian G. fanii differs, however, by having an inversion, trnS-trnR-trnY (TYPE2).
The pairwise divergences in mitogenome data ranged from 0.07 to 3.4% within G. fanii, which is higher than  Length, variable site, nucleotide diversity, and indel. V, variable site; S, singleton; PI, parsimony informative site, π nucleotide diversity; SD, standard deviation; h, number of haplotype; a aligned length; *excluded CNU066493 due to the rearrangement of tRNAs.
that in Pterocladiella capillacea (S.G.Gmelin) Santelices and Hommersand from the Galápagos Islands and Pacific Mexico (Boo et al., 2016a). However, pairwise divergence in rbcL was 0-0.5%, a value lower than interspecific divergences (0.8-0.9%) in the Gelidiales, with the exception of a closely related sibling species . The pairwise divergences in rbcL fell within the large range of reported rhodomelacean divergences, e.g., 0.4-1% in Pairwise divergence within or between species is a relative value rather than an absolute criterion for species recognition. The Australian G. fanii, differed by 2.4-2.8% in cox1 and 0.3-0.5% in rbcL from other populations of G. fanii, may represent a closely related sibling species. This hypothesis may be supported by the unique tRNA inversion. However, because the Long Reef specimen was used in its entirety for the molecular analyses, additional material was not available for further morphological examination to determine whether it is a closely related cryptic species or not (Huisman et al., 2018).
Phylogenomic analyses corroborate the monophyly of G. fanii and its independence from outgroup species. The close relationship of G. fanii with G. acerosa and G. flabella has been supported by previous multigene phylogenies (Boo et al., 2016b,c). Gelidiella fanii is sympatric with G. acerosa in Juan de Nova Island, Tetiaroa Atoll, and Vietnam (this study) as well as Southeast Asia (Lin and Freshwater, 2008;Wiriyadamrikul et al., 2010). Because asexual propagation via tetrasporangia or fragmentation is apparently the only means of reproduction in the genus Gelidiella, we speculate that G. fanii may have arisen by sympatric asexual speciation from an ancestor of the G. acerosa complex. Alternatively, it is possible that sexual reproduction was lost in the genus after speciation, reinforcing the isolation of populations.
The mitogenome and cox1 sequences demonstrate the widespread distribution of G. fanii from Kenya and Juan de Nova Island (Mozambique Channel) to Hawai'i and Tetiaroa Atoll (this study) as well as in northern Australia and Southeast Asia (Lin and Freshwater, 2008;Wiriyadamrikul et al., 2010;Boo et al., 2016c;Huisman et al., 2018). Its occurrence in Eritrea is confirmed by a short, conserved region of rbcL from an archival herbarium specimen. We found several archival herbarium specimens from Samoa and Tonga housed in UC that were morphologically similar to G. fanii; however, our efforts to isolate genomic data from those specimens were unsuccessful.
The broad range of G. fanii, spanning the east and west margins of the Indo-Pacific Ocean, is part of an interesting biogeographical story. In contrast, recent molecular studies have revealed that many so-called widespread species are introduced, misidentified species, or consisted of a complex of local cryptic species (Won et al., 2009;Sherwood et al., 2010;Boo et al., 2018;Díaz-Tapia et al., 2018;Leliaert et al., 2018). That said, this and previous studies have recently distinguished G. fanii from its very similar and widespread congener, G. acerosa, with which it has been long confused.
The phylogeny and haplotype network of Gelidiella fanii detected strong signals of genetic differentiation into seven geographical groups; i) Southeast Asia, ii) Kenya/Juan de Nova Island, iii) Indonesia, iv) northern Australia, v) the Philippines, vi) Tetiaroa Atoll, and vii) Hawai'i. The origin of these groups may be interpreted via two scenarios. One hypothesis is that ancient populations of G. fanii were widely distributed in the Tethys Ocean and the local populations accumulated genetic changes at the geographical periphery. Subsequently, local populations have likely been geographically isolated by changes in sea level and other oceanographic conditions during the Pleistocene epoch (2.6 Ma -11.7 ka), which experienced many glacial-interglacial cycles (Herzschuh et al., 2016). A similar hypothesis has been proposed for Gelidiophycus G.H. Boo, J.K.Park and S.M.Boo and Mazzaella G.De Toni (Montecinos et al., 2012;. This hypothesis is supported by low dispersal capacity of G. fanii, a small species lacking flotation and with asexual reproduction. Our second hypothesis is that Southeast Asia, rich in haplotypes, was the center of origin of G. fanii. Populations may have subsequently migrated to marginal regions of the Indo-Pacific Ocean. This hypothesis has been proposed for Portieria Zanardini species (Gigartinales) that originated in the Central Indo-Pacific and then dispersed in all directions, although the origins of species in remote islands were uncertain (Leliaert et al., 2018). Additional sampling in more regions is necessary to differentiate the two hypotheses, but the logistics of sampling are currently difficult. It is unlikely that G. fanii has recently spread via introductions because of the paucity of shared haplotypes among the seven groups. The many missing haplotypes also indicate deep history.
The mitochondrial genome is a significant source of molecular markers at different levels of resolution. The cox1 gene, with the lowest Ka/Ks ratio (0.0189), is considered a good marker at the population level and has a high degree of resolution at the species level (Boo et al., 2016b). The large size (1,602 bp) of cox1 also provides evolutionarily important information with many singletons (58) and haplotypes (8). Its utility as a population marker has been demonstrated for Gelidiophycus species . COI-5P (∼664 bp), a conserved part of cox1 gene, is the well-known DNA barcoding marker for identifying red algal species and populations (e.g., Saunders, 2005;Sherwood et al., 2010;Hu et al., 2015;Freshwater et al., 2017;Leliaert et al., 2018). The haplotype number (7) and nucleotide diversity (0.0175 ± 0.0026) of COI-5P are similar to those of cox1 gene, supporting its utility in population studies. Both cox2 and cox3 can also be useful markers; cox2 (783 bp) contained 31 singletons and six haplotypes, and cox3 (819 bp) contained 36 singletons and  The hierarchical level "group" indicates 7 groups based on the phylogeny and haplotype network of mitochondrial cox1. d.f. = degrees of freedom, ***P < 0.001. seven haplotypes. The cytochrome b gene (cob, 1,146 bp), with 47 singletons and seven haplotypes, has been used as a DNA barcode at the species and population levels (Saunders and Moore, 2013;Yoon et al., 2014). Five intergenic spacers described in this study are candidates for markers to distinguish populations. The cox2-cox3 (158 bp) contained six haplotypes with a nucleotide diversity of 0.0204 ± 0.0036. This spacer region (∼350-400 bp including flanking parts) has been used to identify red algal species and populations (e.g., Zuccarello et al., 1999;Lim et al., 2013;Kamiya and West, 2014;Pezzolesi et al., 2019). Rps11-nad3 (145 bp) may be a suitable marker with its six haplotypes and nucleotide diversity of 0.0235 ± 0.0037. Nad4-nad5, trnW-trnA, and trnMrns contained seven haplotypes; however, nucleotide diversity was higher in trnW-trnA (0.0422 ± 0.0050) and trnM-rns (0.0483 ± 0.0046) than nad4-nad5 (0.0131 ± 0.0024). Suitable primers for amplifying the above spacer regions have yet to be designed and tested for their practical use in population studies.

CONSERVATION AND CONCLUSION
This study, a detailed characterization of the phylogeography of G. fanii using data from the mitogenomes of eight populations, reveals significant genetic structure correlated with geographic distribution. The continued persistence of G. fanii is unclear, considering its rarity and close association with G. acerosa, with which it likely competes for substrate, nutrients and light. There is a critical shortage of gelidioid biomass for the production of high-grade bacteriological and pharmaceutical agar (Callaway, 2015;Santos and Melo, 2018). The harvest of Gelidiella is possible due to the higher productivity of populations in tropical waters compared to those of temperate species of Gelidium (Ganzon-Fortes, 1994;Santos and Melo, 2018). Commercial harvest of G. acerosa for agar impacts G. fanii and, because it is rare, population reduction over time is probable. Minimizing the loss of genetic diversity is a key factor for conserving economically important marine species. Populations with distinct haplotypes can be recognized and managed as independent conservation units. Populations could also be protected within marine reserves or marine protected areas. Systematic management of natural populations can improve biomass yields for the agar industry. Ecological and phenological studies will be needed in the future to detect changes in population sizes and shifts in geographical ranges due to overharvesting and anticipated oceanographic climate change.
Mitogenome data are currently unavailable for population studies of florideophycean red algae, a large and diverse group comprising 7,000 species. Our results represent a new path toward a mitogenomic picture of the evolution and biogeography of widespread red algal species. This study is the first to describe and analyze complete mitogenomes from the family Gelidiellaceae and to test their utility for population-level studies. We have extended the range of G. fanii from Southeast Asia and northern Australia to the eastern and western margins of the Indo-Pacific Ocean, and have suggested that the observed patterns of genetic groups may have resulted from sea level fluctuations or other environmental changes that occurred during the Pleistocene. On the basis of our mitogenome sequences, we propose that intergenic spacer regions as well as PCGs can be explored as suitable markers at the population level for red algae. Markers derived from mitogenomic analyses are important tools for tracking species ranges and developing future conservation strategies for ecologically and economically important red algal species.

DATA AVAILABILITY STATEMENT
The datasets generated for this study can be found in GenBank under the accession numbers, MT742595-MT742604 for mitogenomes, MT758387-MT758400 for mitochondrial cox1, and MT758377-MT758386 for plastid rbcL sequences.

AUTHOR CONTRIBUTIONS
GHB, MTF, SMB, and KAM conceived and designed the study. GHB, MZ, ARS, SMB, and KAM performed the collections and provided samples. GHB and JRH performed laboratory work and the data analyses. All authors participated in the interpretation of the data and writing the manuscript.