Luminescent characteristics and mitochondrial COI barcodes of nine cohabitated Taiwanese fireflies

Specimen collection and habitat

366 specimens of eight adult and one larval species were randomly collected in flight or from vegetation using hand dip nets from two habitats in the suburbs of Taipei, Taiwan—Nankang (25°01′40.4″N 121°38′02.6″E) and Miaoli County, Nanzhuang (24°37′53.5″N 121°01′37.0″E)—at 18:30–19:30 from April to August 2016–2019 (Table 1). After the bioluminescent spectrum/intensity measurement, the specimens were deposited in a laboratory freezer. For DNA extraction, several specimens were then stored at −80 °C. The remaining specimens were stored at −20 °C for species identification.

The environmental conditions before and after the fireflies began flashing and/or flying were investigated to understand what environmental factors may trigger their nocturnal activity. The temperature, relative humidity, and light intensity (lux) of the firefly habitat were recorded during the period before flashing/flying (18:20∼18:40) and the period after the fireflies began flashing/flying (18:30∼18:50) using HOBO U12-012 data loggers (Onset Computer Corp., Bourne, MA, USA) at 10-sec intervals.

Species identification and morphological measurements

The specimens were collected as previously described in Goh, Lee & Wang (2022). The material collected in this study was identified by LJ Wang on the species level through the use of available references (Ballantyne et al., 2013; Chen, 2003; Jeng, Lai & Yang, 2003; Jeng, Yang & Lai, 2003; Jeng et al., 1998). Five morphological characteristics of the specimens (body length, pronotum length, pronotum width, front wing length, and front wing width) were measured (see Table S3) using a dissecting microscope and photographed with a digital video camera as previously described in Goh, Lee & Wang (2022). During the survey, the specimens were chilled on ice. All surveys were completed within two days of the collection. One to five identified specimens were sacrificed and stored at −20 °C in the Biodiversity Research Center, Academia Sinica, Taipei, Taiwan (contact person: TY Wang, tziyuan@gmail.com).

Bioluminescence spectrum/intensity measurement

The wavelength (λ_max) and luminescent intensity (nW/cm²) of the light flashes produced by the living samples were measured by a USB2000+ spectrometer (Ocean Optics) and a PD300 power meter (Ophir), respectively. All surveys were completed within two days of the collection. The wavelength and luminescent intensity measurements were performed in a dark room by directly attaching the detector of the USB2000+ spectrometer or PD300 to the light organ of a trapped firefly (Fig. 1 modified from Goh, Lee & Wang (2022)). The average wavelength peak and λ_max were obtained from an average of 3–5 measurements in complete darkness at 25 °C with 75% humidity. The luminescent intensity of the flash was obtained by averaging each flash from 3–10 min of recording data with a PD300 power meter. To compare the luminescent intensity data from PD300 and HOBO U12-012 using the same units, all data in the energy unit nW/cm² were converted into lux via the conversion 1 lux = 1.464E−07 W/c m² = 146.41 nW/cm² (at 555 nm).

Statistics

The differences in bioluminescence spectrum among the specimens were determined by the Chi-square test between two species.

DNA barcode sequencing

Crude DNA was extracted from thoracic muscles via the ZR Tissue & Insect DNA MicroPrep™ kit (D6015). Two beetle-specific primers (ClepFolF 5′-ATTCAACCAATCATAAAGATATTGG-3′ and ClepFolR 5′-TAAACTTCTGGA TGTCCAAAAAATCA-3′) were designed based on the comprehensive DNA barcode database of beetles (Hendrich et al., 2015) to amplify a 620-bp segment including the cytochrome oxidase I (COI) gene. Polymerase chain reactions (PCRs) in 50-µL volumes were performed with a dNTP concentration of 200 µM and a primer concentration of 0.3 µM, with 50 ng of genomic DNA, one unit of TaKaRa Taq™ DNA Polymerase, and the buffer supplied by the manufacturer. The PCR was run for 35–40 cycles under the following conditions: denaturation at 95 °C for 30 s, annealing at 50∼55 °C for 40 s, extension at 72 °C for 1 min, and a final extension at 72 °C for 10 min. The product mixture was used as a template for DNA sequencing (Genomics Ltd., Taipei, Taiwan). Haplotype sequences were deposited into GenBank under accession numbers MT534191–MT534201, ON209457 (Table 2).

Molecular phylogeny

The COI sequences of closely-related species and/or species with known λ_max of luminescence were downloaded from GenBank based on previous studies (Adamowicz, 2015; Arnoldi, Neto & Viviani, 2010; Cassata, 2020; Dong et al., 2021; Han et al., 2020; He et al., 2021; Hendrich et al., 2015; Jusoh et al., 2021; Jusoh et al., 2018; Jusoh et al., 2014; Kim et al., 2001; Li, Yang & Fu, 2022; Liu et al., 2017; Liu & Fu, 2020; Martin et al., 2017; Martin et al., 2015; Muraji, Arakaki & Tanizaki, 2012; Oba, Branham & Fukatsu, 2011; Oba et al., 2020; Osozawa et al., 2015; Roslin et al., 2022; Rulik et al., 2017; Stanger-Hall, Lloyd & Hillis, 2007; Usener & Cognato, 2005; Wilcox, 2021; Zhang et al., 2018). Sequences were then aligned using the ClustalX program (Thompson, Gibson & Higgins, 2002) with the default setting in MEGA X (Kumar, Stecher & Tamura, 2016), followed by length trimming due to different amplicons. After trimming, the short-length sequences were removed. At most, three representative sequences were kept for each species to simplify the tree topology. There were a total of 520 positions and 161 nucleotide sequences in the final dataset (Table S1). Neighbor-joining (NJ) (Saitou & Nei, 1987) and maximum-likelihood (ML) trees were constructed using GTR+G+I distances in MEGA X with 500 bootstrap replications (Felsenstein, 1985). The substitution model (parameter) used to calculate GTR+G+I distances (Nei & Kumar, 2000) was selected using Modeltest v3.7 (Posada & Crandall, 1998). The differences in the composition bias among sequences were considered in the evolutionary comparisons (Tamura & Kumar, 2002).

Cohabitated species composition at Nanzhuang and Nankang

From April to August 2016–2019, we collected 366 flying specimens from two firefly habitats (Nanzhuang and Nankang) in northern Taiwan (Table 1 & S2 Fig. 2), and morphologically identified them to the species level (Table S3). These specimens comprised nine different species, including adult males and/or females of Aquatica ficta, Luciola filiformis, Abscondita cerata, Luciola kagiana, Luciola curtithorax, Abscondita chinensis, Curtos sauteri, and Curtos costipennis. Only Pyrocoelia praetexta was observed as larvae from April to August. Five of the species—Aq. ficta, L. filiformis, Abs. cerata, L. kagiana and P. praetexta—were found in both habitats. Luciola curtithorax was collected only in Nankang, and Abs. chinensis, C. sauteri, and C. costipennis were collected only in Nanzhuang (Table S2). To simplify the results, Table 1 combined all specimens from the two habitats for further analysis.

Occurrence periods of the cohabitated species

Figure 3 shows the estimated occurrence periods of the eight adult species in Nankang and Nanzhuang based on the collection dates of the specimens. Abscondita cerata and L. kagiana occurred in April–May, while C. sauteri and Aq. ficta occurred in May–August. L. curtithorax occurred in May-July while L. filiformis occurred in May–June. Abscondita chinensis was found only in June, while C. costipennis was found only in August. Six of the eight studied species shared an overlap occurrence period in May: C. sauteri, Aq. ficta, Abs. cerata, L. kagiana, L. filiformis, and L. curtithorax. The occurrence periods of this study overlapped with the occurrence periods of previous studies based in other habitats in Taiwan (Chen, 2003; Chen & Jeng, 2012).

Differences in the luminescence spectrum between the cohabitated species

The average λ_max from the luminescent spectra of nine studied species (eight adult species and larval P. praetexta) ranged from about 552 nm (green-yellow) to 575 nm (yellow-orange) (Table 1, Fig. S1). Excluding insufficient data on three species, the average λ_max of five species showed no significant difference between intraspecific males and females. The average λ_max of L. curtithorax was significant different between female and male. Ignoring the sexual differences, the pairwise comparison of interspecific λ_max (Table 3) showed that the studied species commonly displayed significant difference (p-values <0.05) in λ_max to most other studied species. No significant difference was found between Abs. chinensis and L. curtithorax (p-value = 0.0604) or between C. sauteri and P. praetexta (p-value = 0.161).

To determine the courtship behaviors of cohabited fireflies based on sex, we further compared the λ_max between six interspecific adult females (Table 4) and males (Table 5). Most studied species revealed significant differences in λ_max between interspecific females. However, no significant difference has found between females of Aq. ficta and Abs. cerata (p-value = 0.0742). In contrast, the λ_max comparison between interspecific adult males (Table 5) showed no significant difference between Aq. ficta and Abs. cerata (p-value = 0.672), between L. curtithorax and Abs. chinensis (p-value = 0.626), between L. kagiana and Abs. chinensis (p-value = 0.1835), and between L. kagiana and L. curtithorax (p-value = 0.2606).

Correlation between firefly luminescent intensity and environmental photic intensity

This study was performed during an Abs. cerata massive occurrence (April to May) in Nankang and Nanzhuang. During the studied periods, it was estimated that the average environmental light intensity during twilight, the ten-minute period before the fireflies started flashing or flying in the habitats, was in a range of 35.7–136.5 lux (Table 6). The suitable environmental light intensity for fireflies flashing and/or flying was in a range of 6.49–28.1 lux.

With the exception of L. kagiana and P. praetexta due to abnormal behavior (no glowing or glowing in extremely low light intensity), the luminescent intensity of seven adult species was about 1.2–14 lux (182.1–2,048 nW/cm²) in male fireflies and nearly 0.8–5.8 lux (122.8–850 nW/cm²) in female fireflies (Table 1). The results showed that the male fireflies have higher luminescent intensity than the females, which might be related to their courtship behaviors.

Herein, we argue that firefly luminescent intensity is correlated with environmental photic intensity. For examples, among the studied male species, male A. cerata produced the brightest flashes, measuring up to 14 lux (or 2,048 nW/cm²). In contrast, female A. ficta emitted the brightest flashes among the studied female species, measuring up to 5.8 lux (or 850 nW/cm²). In addition, the maximum luminescent intensity emitted from the five kinds of adult males was 2.3–14 lux (332.1–2,048 nW/cm²), which is 1.01–7.26-fold higher than that of conspecific females (1.9–5.8 lux or 282.2–850 nW/cm²). Thus, this result clearly shows that the range of the environmental light intensity (6.49–28.1 lux) when fireflies begin to flash partially overlaps with the luminescent intensity of fireflies. In addition, during 18:00–19:30, the change in average environmental temperature and relative humidity were in a range of 17.1–25.0 °C and 71.2–95.8%, respectively.

Molecular phylogeny of Lampyridae inferred by COI barcodes

To reveal how bioluminescence evolved, it is important to compare the luminescence spectrum and molecular phylogeny. The COI barcodes of eight studied species (except C. costipennis) were successfully sequenced for phylogenetic analysis (Table 2). All are new COI barcodes of Taiwanese fireflies sequenced in this study. Their haplotype sequences were deposited in GenBank under accession numbers MT534191 –MT534201, ON209457.

The NJ tree (Fig. S2) and ML tree (Fig. S3) indicate that the studied genera Abscondita, Curtos, Aquatica, and Luciola belong to Luciolinae, while the genus Pyrocoelia belongs to Lampyrinae, a monophyly supported by previous mitogenomic phylogeny (Wang, Wu & Wang, 2021). However, the short COI sequences showed incongruence grouping among subfamilies in the high-level phylogeny. For example, Rhagophthalmus (Rhagophthalmidae) was placed close to the Luciolinae with a low bootstrap value; Stenocladius did not form a clade with Drilaster as Ototretinae. Nevertheless, most studied species are placed correctly with congeners (Figs. 4 and 5).

There are several monophyletic clades supported by medium or high bootstrap values. Lampyrinae was a monophyletic clade with Pyrocoelia, Diaphanes, Lampyris, Microphotus, Aspisoma, Photinus, Pleotomodes, Ellychnia, Lucidina, Phosphaenus, Pyractomena, and Pyropyga, although these genera did not form stable sister groups with each other. Photuris and Bicellonycha (as Photurinae) formed a clade with Pyropyga nigricans.

A previous study indicated that Luciolinae is not monophyletic, even when 436 gene loci were used (Martin et al., 2019). Thus, it is reasonable to see polyphyletic Luciolinae in the COI gene tree. The Luciolinae complex included the monophyletic genera Luciola, Aquatica, Pterophyx, Sclerotia, Abscondita, Pygoluciola, Curtos, and Lamprigera, which comprises Pristolycus, Vesta, and Emeia pseudosauteri. There are still several monophyletic clades with medium bootstrapping values. Excluding L. cruciata, the COI barcode grouped 12 Luciola species as a monophyly supported by a medium bootstrapping value (61/54), including the type species (L. italica). Luciola cruciata formed a stable clade with five Aquatica species. In addition, Pygoluciola clustered with Abscondita while Curtos clustered with Lampyrigera supported only with a low bootstrapping value.

Large COI sequence variations can be found between and/or within geographically distinct species. For example, two COI sequences of Chinese Emeia pseudosauteri were separated into distinct clades. Emeia pseudosauteri is restricted to central China and isolated among mountains. Such habitat isolation caused great mitochondrial DNA variation (Liu & Fu, 2020). Accordingly, there might be cryptic species and a need to reclassify some other species. A detailed analysis of Taiwanese fireflies will be discussed later.

Differences in luminescence spectra among cohabitated fireflies

The λ_max values of the luminescence spectra were similar within the same species, but different between species (Table 3). The pairwise comparison also showed significantly different interspecific λ_max, except between adult Abs. chinensis and adult L. curtithorax and between adult C. sauteri and larval P. praetexta. A recent study (Goh, Lee & Wang, 2022) recorded species-specific flash patterns of three sympatric male fireflies (Abs. cerata, L. kagiana, and L. curtithorax); both that study and the present study (Table 5) showed that males of L. kagiana and L. curtithorax have similar λ_max values but still retain their own unique flash patterns. Thus, the four studied species with similar λ_max might also have species-specific flash patterns, which need further study in the future. The above results might imply that most cohabitating fireflies distinguish between each other based on different luminescence spectra and/or specific flash patterns. This implication is important to consider when previous literature (Chen, 2003; Chen & Jeng, 2012) and our findings indicate that all studied adult fireflies appear simultaneously from April to June (Fig. 3). Thus, various cohabitated species may have evolved species-specific recognition to improve male–female searching within such a densely populated area over such a short nocturnal activity time.

Most fireflies have significantly different λ_max between interspecific females (Table 4). Only those of Aq. ficta and terrestrial Abs. cerata females were not significantly different. However, the microhabitat of Aq. ficta and Abs. cerata was in the aquatic habitat and moist forest, respectively (Chen, 2003; Jeng, Yang & Lai, 2003). The flying males could still have better chance to find the conspecific females in their specific microhabitat.

In contrast, similar λ_max were found in males of L. curtithorax, Abs. chinensis, and L. kagiana (Table 5), but males of L. curtithorax and L. kagiana have their own flash patterns (Goh, Lee & Wang, 2022). Thus, flash pattern is another key for cohabited female fireflies to recognize conspecific males (Lewis & Cratsley, 2008; Lower, Stanger-Hall & Hall, 2018). In addition, Abs. cerata and L. kagiana have different nocturnal activity time, while L. curtithorax is restricted to the dark ground layer of forest (Goh, Lee & Wang, 2022).

Based on the above phenomena, different λ_max, species-specific flash patterns, microhabitat choices, nocturnal activity time, and/or isolated mating seasons are key factors that may lead to the species-specific courtship of cohabitated fireflies.

Luminescent intensity of flashes implies sensing distance

Fireflies seem to be very sensitive to the photic environment in the evening. Artificial light pollution is a major force influencing firefly proliferation, mating, and growth (Costin & Boulton, 2016; Firebaugh & Haynes, 2016; Haynes & Firebaugh, 2019; Owens, Meyer-Rochow & Yang, 2018). The environmental light intensity and the light sensitivity of the fireflies influence whether the fireflies will flash. Therefore, the luminescent intensity of flashes emitted by fireflies could be an ecological indicator for evaluating light pollution to fireflies. This study further investigated this issue based on the first flash time of abundant Abs. cerata (Table 6). Fireflies start flashing or flying (nocturnal activity) when the environmental light intensity decreases to 6.49–28.1 lux (∼950–4,114 nW/cm²). The luminescent intensity of male Abs. cerata ranges from the average (2.1–3.4 lux or 406.6 ± 96.5 nW/cm²) to the maximum (14 lux or 2,048 nW/cm²), which overlaps with the environmental light intensity suitable for their nocturnal activity (Table 1). A previous study also showed that the abundant Abs. cerata begins flashing when the photic environment decreases to 0.04–1.38 lux (Ohba & Yang, 2003). All imply that Abs. cerata could tolerate environmental light intensity around 28.1 lux but wait until 6.49 lux to start nocturnal activity in the evening at twilight. Another study also revealed that most male Abs. cerata start to fly in the evening at twilight while L. kagiana starts its nocturnal activity later (Goh, Lee & Wang, 2022), which Table 1 indeed showed lower luminescent intensity of L. kagiana. In addition, another study (Owens, Meyer-Rochow & Yang, 2018) revealed that half of the Aq. ficta specimens stopped flashing under bright exposure (∼20 and 200 lux). Table 1 further shows that the luminescent intensity of male Aq. ficta ranged from the average (3.1–4.1 lux or 525.7 ± 71.1 nW/cm²) to the maximum (7.5 lux or 1102 nW/cm²), which we also observed the Aq. ficta appeared with Abs. cerata during the same period of nocturnal activity in Nanzhuang. Such differences in luminescent intensity of the three species might imply another adaptation factor for the different nocturnal activity time among species.

Next, we measured the putative sensing distance between males and females. During a typical courtship, the flying males flash to attract perched females. Then, the female responds and flashes to the flying male. The male fireflies close and lands near the female; each displays different flash patterns for communication. As they court each other, the paired fireflies stop flashing on perch. Communication between female and male fireflies relies on the illumination of their light organ in the dark. Usually, the average luminescent intensity emitted by most females (the light organ from single tagma) is around half that of males (the light organ from double tagmata). The differences in luminescent intensity between sex could be due to their courtship behavior for sensing each other. The male needs a higher intensity exposure for females to find him while the female needs to save energy for later proliferation and only responds to male signals with detectable intensity.

The sensing distance between a female and male could be relative to their bioluminescent intensity. So, using the luminescent intensities of male and female, we could estimate the sensing distance. The assumption is the females have higher sensitivity while males have higher luminescent intensity. So, the luminescent intensity difference between male and female could be the sensing ability for a female to detect a male or vice versa. Thus, the maximum luminescent intensity might represent the maximum sensing distance between females and males, assuming that the minimum sensing distance (r, meter) is around the same luminescent intensity between females and males.

We can estimate the sensing distance using the example of the Abs. cerata. The males have a maximum luminescent intensity of 2,048 nW/cm² (14 lux) and the females have a maximum luminescent intensity of 282.2 nW/cm² (∼1.93 lux). Using the formula 14 / (r²) = 1.93, we can estimate the maximum sensing distance (r) for this species to be around 2.7 m. Using the same formula calculation with average luminescent intensities, we estimated the average sensing distance to be around 1.8 m. In other words, the putative sensing distance for female Abs. cerata could range from 1.8 to 2.7 m, which may also be the sensing distance for a flying male searching for a female. That said, it is important to note that most females prefer to perch as males fly to approach them (Goh, Lee & Wang, 2022), since perched females should flash less than what we measured. Thus, the sensing distance between males and females may actually be shorter. Nevertheless, the luminescent intensity could be an indicator of the sensing distance between flying males and perched females. After all, a previous study revealed that male Photinus carolinus use a 15–30 cm landing distance when approaching perched females (Copeland, Moiseff & Faust, 2008), which is a reasonable sensing distance in our estimation. Further behavior experiments should investigate these issues.

Monophyly of Luciola sensu stricto

Both the mitogenome (Jusoh et al., 2021) and COI barcode (this study) revealed that each of the studied L. species form a clade, except for L. cruciata. Luciola cruciata and genus Aquatica were grouped together. The other genera of Luciolinae (Curtos, Pteroptyx, Sclerotia, Abscondita, Pygoluciola) are distinct.

Lamprigera is not within Lampyrinae

Both the mitogenome (Wang, Wu & Wang, 2021) and 436 nuclear loci (Martin et al., 2019) indicated genus Lamprigera groups within Luciolinae instead of Lampyrinae. The COI phylogeny (Figs. 4 and 5) also showed that Lamprigera is a sister group to Curtos and separate from Lampyrinae. In addition, the morphology and COI sequences of eight native species (Dong et al., 2021) further revealed that Lamprigera should be closer to Luciolinae.

Cryptic species implied by mitochondrial COI barcode variation

The mitochondrial genetic variation of fireflies within a population or adjacent regions has fewer genetic differences—e.g., the desert-based Microphotus octarthrus (Usener & Cognato, 2005), the widespread Photinus pyralis (Lower, Stanger-Hall & Hall, 2018), and the Korean Aquatica lateralis (Kim et al., 2001; Suzuki et al., 2004). Previous biogeographical study revealed that the two studied sites in Northern Taiwan are within the same geographical regions; thus, we sequenced 1–3 individuals, except the abundant Abs. cerata, in which only one SNP site between two haplotypes could be found from nine Abs. cerata individuals of two habitats. Herein, the COI barcode showed a genus-level resolution for species identification in Figs. 4 and 5, although COI phylogenies in higher-level topologies are not consistent with those of previous morphological studies (Ballantyne et al., 2013; Ballantyne et al., 2015; Ballantyne & Lambkin, 2013; Martin et al., 2017; Stanger-Hall, Lloyd & Hillis, 2007) and molecular phylogeny (Chen et al., 2019; Martin et al., 2017; Martin et al., 2019; Wang, Wu & Wang, 2021). Nevertheless, the COI barcode could successfully identify most species at the genera-to-species level (Figs. 4 and 5). The COI phylogeny showed that the studied genera Abscondita, Curtos, Aquatica, and Luciola belong to Luciolinae, while Pyrocoelia belongs to Lampyrinae as expected.

The COI sequence variations revealed several cryptic species in Taiwan. For example, 62 SNP sites (∼11.9% variation) were found in the COI sequences between Taiwanese and northern Chinese Abs. chinensis. Sixty-four SNP sites (∼12.3% variation) were found in the COI sequences between Abs. terminalis and Chinese (Taiwanese) Abs. chinensis, respectively. Building on a previous study (Ballantyne et al., 2013), this study further showed that the λ_max of the Abs. chinensis lantern spectrum (flash color) is different between Taiwanese (572 nm) and northern Chinese (565 nm) individuals, though there might be unknown environmental effects that cause the flash color variation in widespread species, like with the North American firefly, Photinus pyralis (Lower, Stanger-Hall & Hall, 2018). Thus, the new evidence reveals that Taiwanese Abs. chinensis may be a distinct species to Chinese Abs. chinensis (Figs. 4D and 5D).

Large COI variation was also found in six Asian species: P. praetexta (Figs. 4A and 5A), C. costipennis (Figs. 4E and 5E), Aq. ficta and Aq. lateralis (Figs. 4C and 5C), and L. curtithorax and L. filiformis (Figs. 4B and 5B). The COI barcode also indicated 17 SNP sites (∼3.27% variation) between Taiwanese and northern Chinese Aq. ficta (Figs. 4C and 5C). One study indicated that the characterization of the Chinese Aq. ficta differed slightly from the Taiwanese Aq. ficta (Ballantyne & Lambkin, 2009). In addition, Aq. leii was considered as a different species to the Chinese Aq. ficta (Fu, Ballantyne & Lambkin, 2010). However, there are only two SNP sites (∼0.39% variation) between Chinese Aq. ficta and Aq. leii. The Taiwanese Aq. ficta may be a cryptic species with a large variation (17 SNP sites). In contrast, the COI barcode indicated a large variation of 10% (52 SNP sites) between Korean and Japanese Aq. lateralis (Figs. 4C and 5C). There are 58 SNP sites (∼11.2% variation) between Taiwanese and northern Chinese L. curtithorax (Figs. 4B and 5B). The COI barcode indicated 58 SNP sites (∼11.2% variation) between Taiwanese L. filiformis and Japanese L. filiformis yayeyamana (Figs. 4B and 5B). The COI barcode indicated 57 SNP sites (∼11% variation) between Taiwanese and southwestern Chinese P. praetexta (Figs. 4A and 5A). The COI barcode indicated 33 SNP sites (∼6.35% variation) between southern Japan and eastern Chinese C. costipennis (AB608764 and MK609965 in Figs. 4E and 5E). All these examples indicate large COI variations between two geographical isolates. Further investigations are needed to reclassify these geographically isolated species.

Bioluminescent evolution inferred from mitochondrial COI barcodes and known phylogeny of Lampyridae

The contracted high-level phylogeny (Chen et al., 2019; Martin et al., 2019) and Luciolinae grouping (Jusoh et al., 2021) correspond well with our bioluminescent evolution phylogeny (Fig. 6). Another study revealed the bioluminescent evolution via recombinant luciferases and suggested the origin of beetle bioluminescence (Oba et al., 2020). Accordingly, this study gives a detailed summary on the evolution of bioluminescence in Lampyridae based on the λ_max of its luminescence spectrum (Arnoldi, Neto & Viviani, 2010; Goh, Lee & Wang, 2022; He et al., 2021; Oba et al., 2020; Wilcox, 2021). Our studied species further revealed that the fireflies’ luminescence color was originally a green color in a Lampyridae ancestor, then red-shifted to a yellow-green in Luciolinae and is now an orange-yellow color in some derived species (Fig. 6).