Genetic Architecture of the Variation in Male-Specific Ossified Processes on the Anal Fins of Japanese Medaka

Traits involved in reproduction evolve rapidly and show great diversity among closely related species. However, the genetic mechanisms that underlie the diversification of courtship traits are mostly unknown. Japanese medaka fishes (Oryzias latipes) use anal fins to attract females and to grasp females during courtship; the males have longer anal fins with male-specific ossified papillary processes on the fin rays. However, anal fin morphology varies between populations: the southern populations tend to have longer anal fins and more processes than the northern populations. In the present study, we conducted quantitative trait locus (QTL) mapping to investigate the genetic architecture underlying the variation in the number of papillary processes of Japanese medaka fish and compared the QTL with previously identified QTL controlling anal fin length. First, we found that only a few QTL were shared between anal fin length and papillary process number. Second, we found that the numbers of papillary processes on different fin rays often were controlled by different QTL. Finally, we produced another independent cross and found that some QTL were repeatable between the two crosses, whereas others were specific to only one cross. These results suggest that variation in the number of papillary processes is polygenic and controlled by QTL that are distinct from those controlling anal fin length. Thus, different courtship traits in Japanese medaka share a small number of QTL and have the potential for independent evolution.

ABSTRACT Traits involved in reproduction evolve rapidly and show great diversity among closely related species. However, the genetic mechanisms that underlie the diversification of courtship traits are mostly unknown. Japanese medaka fishes (Oryzias latipes) use anal fins to attract females and to grasp females during courtship; the males have longer anal fins with male-specific ossified papillary processes on the fin rays. However, anal fin morphology varies between populations: the southern populations tend to have longer anal fins and more processes than the northern populations. In the present study, we conducted quantitative trait locus (QTL) mapping to investigate the genetic architecture underlying the variation in the number of papillary processes of Japanese medaka fish and compared the QTL with previously identified QTL controlling anal fin length. First, we found that only a few QTL were shared between anal fin length and papillary process number. Second, we found that the numbers of papillary processes on different fin rays often were controlled by different QTL. Finally, we produced another independent cross and found that some QTL were repeatable between the two crosses, whereas others were specific to only one cross. These results suggest that variation in the number of papillary processes is polygenic and controlled by QTL that are distinct from those controlling anal fin length. Thus, different courtship traits in Japanese medaka share a small number of QTL and have the potential for independent evolution. In the males of several fishes, particular rays of the median fins are modified for attracting females and/or copulation during courtship (Meyer 1997;Zauner et al. 2003;Kwan et al. 2013). For a particular body part to be modified, that part should be genetically and developmentally decoupled from other body parts because the presence of strong genetic correlations would lead to correlated evolution of other parts (Lande 1979), which is potentially deleterious (Badyaev 2004). Therefore, investigating the genetic architectures underlying variation in courtship traits is indispensable to obtain a better understanding of the evolutionary mechanisms underlying the diversification of courtship traits. Genetic mapping of courtship traits has been extensively conducted in insects, and these studies have demonstrated that quantitative trait loci (QTL) controlling different components of courtship traits occasionally, but not always, overlap (Howard and Medaka fishes of genus Oryzias provide a good system for investigations into the genetic architecture underlying courtship traits because sexually dimorphic traits, including modified median fins, exhibit variation between closely related species (Iwamatsu 2006). In Oryzias, males generally have longer anal and dorsal fins than females, although the magnitude of sexual dimorphism differs between species and even populations (Iwamatsu 2006;Kawajiri et al. 2009Kawajiri et al. , 2014. For example, southern Japanese medaka (Oryzias latipes) males have longer anal and dorsal fins than northern males (Kawajiri et al. 2009;Kawajiri and Yamahira 2011). Experiments examining Japanese medaka indicate that longer anal fins are likely to increase reproductive success in males. During courtship of Japanese medaka, the male and the female are aligned side-by-side in the same direction, and the male wraps the female with his dorsal and anal fins; then, the female spawns her eggs and the male ejaculates his sperm. Previous research indicated that experimental shortening of male anal fins reduced the fertilization rate (Egami and Nambu 1961;Koseki et al. 2000). Additionally, a behavioral study showed that females preferred males with longer fins .
Males of several medaka species also possess male-specific ossified processes called papillary processes ( Figure 1); in Japanese medaka, the expression of these processes is induced by androgens (Ogino et al. 2013;Nakamura et al. 2014). The presence and the number of papillary processes differ between species and populations. For example, all Japanese medaka males possess papillary processes; however, medaka species from Indonesia lack papillary processes, although Indonesian males have longer anal fins than females. Papillary processes are thought to rub the female to induce spawning (Ono and Uematsu 1957), although experimental tests of the functions of papillary processes have not yet been conducted. Furthermore, Japanese medaka females sometimes reject the release of eggs by disentangling themselves from the male's anal fins (Fujimoto et al. 2015), and the papillary processes in the males may prevent the females from escaping. Thus, the length of the anal fin and the presence of papillary processes show different patterns of variation in the genus Oryzias.
The present study had two main aims. First, we conducted QTL mapping of the number of papillary processes on the anal fin in an F 2 family that was used previously for QTL mapping of the anal fin length (F 2 intercross derived from a southern Japanese medaka female and a northern Japanese medaka male) . This result was used to test whether the QTL for the papillary process numbers overlapped with the previously identified QTL for the anal fin length. Second, we made another independent cross derived from a northern Japanese medaka female and a southern Japanese medaka male and conducted QTL mapping of both the papillary process numbers and the anal fin length. This result was used to test the number of QTL that were replicated between two independent F 2 crosses. Because multiple alleles with different effects are likely segregating within natural populations (Mauricio 2001), QTL mapping with multiple crosses is necessary to understand the genetic variation underlying phenotypic variation in natural populations (Barnwell and Noor 2008).

Fish collection and crosses
Pictures of wild-caught adult medaka fishes that were taken previously were used for the analysis of wild-caught fishes. The southern population was collected from Ginama, Okinawa Prefecture, Japan (26°839N, 128°2 69E) in 2012, whereas the northern population was collected from Mayajiri, Aomori Prefecture, Japan (40°509N, 140°499E) in 2006 and 2011. In the laboratory, 22 southern males and 20 northern males were photographed individually from the side with the use of a Pentax Optio W90 (Ricoh Imaging, Tokyo, Japan) and Coolpix 4500 (Nikon, Tokyo, Japan), respectively. The papillary processes on each anal fin ray of each male were counted with the photographs.
To make F 2 crosses, the same populations were collected in 2011. One F 2 cross was performed with an Okinawa female and an Aomori male (termed OFAM) that had been used for QTL mapping of the anal fin length, as reported elsewhere . This family consisted of 78 males and 92 females. We performed another F 2 cross with an Aomori female and an Okinawa male (termed AFOM). This new cross was performed according to previously described procedures . To summarize, one female from the north and one male from the south were crossed to create F 1 progeny. Next, F 1 females and F 1 males were allowed to mate to produce an F 2 family composed of 102 males and 87 females. The F 2 progeny were grown from the egg stage individually in single acrylic containers (50 · 100 mm with a depth of 200 mm) and maintained under a 14:10-hr lightdark cycle.

Morphological analysis of the F 2 crosses
Each fish was photographed at 11 ages: 15,20,27,34,43,52,63,74,89,104, and 124 d after fertilization (DAF). For photographing, the fish were placed in a narrow transparent acrylic box and photographed from the side with a digital camera (Pentax Optio W90, Ricoh Imaging, Tokyo, Japan) as described previously (Kawajiri et al. 2009). The papillary processes on each anal fin ray were counted with only the photographs taken at DAF 124, because fish at DAF 104 and younger were too small to allow counting of the papillary processes based on the pictures.
The standard length and anal fin length of each fish were measured at each age to the nearest 0.01 mm with Adobe Illustrator 10.0.3 (Adobe Systems, San Jose, CA) with the plug-in software BPT-Pro2 (Baby Universe, Fujisawa, Kanagawa, Japan). The anal fin length was measured at the fourth anal fin ray, as described previously (Kawajiri et al. 2009. In addition to the fin length at each age, we measured the growth curve shapes of the anal fins using orthogonal polynomial regression. Details regarding the methods were described previously , and the R scripts used are available from Dryad (http://dx.doi.org/10.5061/dryad.nq36d). In brief, a cubic spline was fitted to the mean of the growth trajectories of all F 2 fish, and the residuals from the mean were calculated for each F 2 individual. Third-order Legendre polynomials were fitted to these residuals. Then, we derived the weights of orthogonal polynomials from each fit to describe the curve shape of each F 2 fish; these weights (from order 0 through order 3) were used for QTL mapping.

Linkage map construction and QTL mapping
The linkage map of the OFAM family was previously reported . The linkage map of the AFOM family was constructed according to previously described procedures . In brief, genomic DNA was isolated from the right pectoral fin and the right eye of each F 2 fish with a DNeasy Blood &Tissue Kit (QIAGEN, Valencia, CA). Single-nucleotide polymorphism (SNP) genotyping of the grandparents and F 2 progeny was conducted with the custom VeraCode GoldenGate Genotyping Assay with 384 SNP panels (Illumina, San Diego, CA) that were previously designed for Japanese medaka fish . Details regarding the design of the SNP panel are available from Dryad (http://dx.doi.org/10.5061/dryad. nq36d). The data were analyzed with Illumina GenomeStudio software, and 273 successfully genotyped SNPs were used for linkage map construction. JoinMap 3.0 software was used to construct a linkage map with a LOD threshold of 4.0 (Van Ooijen and Voorrips 2001).
QTL mapping of the papillary process number was conducted by multiple QTL mapping in R/qtl (Broman and Sen 2009;Arends et al. 2010). Unsupervised backward elimination was used to find cofactors. Figure 2 Boxplots of the number of papillary processes. The total number of papillary processes (Total) and the number of papillary process on each fin ray (Ray) were compared between the southern males (S) and the northern males (N). The asterisks indicate significant differences in the traits between crosses by the Mann-Whitney U-test (P , 0.001). N = 22 for the southern population and N = 20 for the northern population. Because only two southern males had papillary processes on Ray 10 and only one northern male had a papillary process on Ray 18, these two fin rays were not included in the analysis.
n Genome-wide permutation tests (1000 bootstraps) were conducted to find significant QTL (P , 0.05). Significant QTL were used to construct a final model, and the calculation of the percentage of variance explained by each QTL and F-tests of the effects of QTL were conducted with the fitqtl function in R/qtl. Because multiple QTL mapping in R/qtl cannot be used to analyze X chromosomes, we excluded the X chromosome (LG1) from the analysis. Instead, we conducted a genome scan with a single QTL model (the function scanone in R/qtl), but we did not find any significant QTL on LG1. Therefore, we determined that the exclusion of LG1 did not affect our conclusions. The dominance of each QTL was analyzed with the function effectscan in R/qtl. A |d/a| ratio . 0.5 was considered dominant, while a |d/a| ratio , 0.5 was considered additive (no QTL had a |d/a| ratio . 1.0) (Bernatchez et al. 2010). The genotype and phenotype data that can be used to reproduce the results are available from Dryad (doi: 10.5061/dryad.5v71p).

Statistical analysis
A Mann-Whitney U-test was conducted to compare the number of papillary processes between populations using the exactRankTests R package (Hothorn and Hornik 2015). Pearson's correlation coefficients were calculated to test the correlation between traits. P-values of LOD scores from the QTL mapping were determined with 1000 genomewide bootstraps. When the 95% Bayesian credible intervals overlapped, two QTL were considered to overlap.

Variation in the number of papillary processes
We first compared the number of papillary processes of wild-caught fish between populations (Figure 2). The southern males had a significantly greater total number of papillary processes than the northern males (Mann-Whitney U-test, W = 414, P , 0.001). Then, we examined the number of papillary processes on each fin ray separately. We found that the number of papillary processes on Rays 12216 significantly differed between populations, with the southern males having more processes than the northern males (Mann-Whitney U-test, P , 0.001).
Genetic mapping of the variation in the number of papillary processes in the OFAM family QTL mapping of the total number of papillary processes in the OFAM family revealed two significant QTL (Table 1, Figure 3, and Supporting Information, Figure S1). A comparison of phenotypic traits among three genotypes (i.e., southern homozygotes, heterozygotes, and northern homozygotes) at these loci showed that the QTL on LG11 had an expected direction of effect such that the southern allele increased the total number of papillary processes (Figure 4). In contrast, the QTL on LG17 had the opposite direction of effect: the southern allele reduced the total number of papillary processes (Figure 4). The southern allele was dominant at the QTL on LG11 (|d/a| = 0.892), whereas the northern allele was dominant at the QTL on LG17 (Figure 4, Table 1; (|d/a| = 0.822)). Next, we conducted QTL mapping of the number of papillary processes on each fin ray separately. We found at least one significant QTL for the number of papillary processes on Rays 13217, whereas no significant QTL were found for Ray 11 or Ray 12 (Table 1, Figure 3, and Figure S1). In addition to LG11 and LG17, we found significant QTL on LG14, LG19, and LG22, suggesting that these QTL were specific to particular fin rays and contribute little to the total number of papillary processes. A comparison of phenotypic traits among genotypes of these fin ray2specific QTL showed that all of these QTL on LG14, LG19, and LG22 had the expected directions of effects ( Figure 4). However, the degree of dominance varied among QTL with |d/a| ranging from 0.152 to 0.953 (Table 1). A significant interaction was only found between LG11 (24 cM) and LG19 (34 cM) for Ray 15 (P-value of F-test of the interaction = 0.0249, percentage of variance explained by interaction = 9.6): the QTL on LG11 became effective in the presence of the northern allele at LG19 ( Figure S2).
A heat map of the LOD scores revealed that different fin rays were controlled by different QTL, although nearby fin rays tended to share QTL ( Figure 5). For example, the profiles of the LOD scores appeared to be similar either between Ray 12 and Ray 13, between Ray 14 and Ray 15, or between Ray 16 and Ray 17. However, the profiles of the LOD scores of Rays 12-13 appeared to differ compared to those of Rays 16217 (see also Figure S1). These results indicate that the numbers of papillary processes on different parts of the anal fins are controlled by different QTL.
Because the number of papillary processes was correlated with the standard length and anal fin length (Table S1), we conducted QTL mapping with either the standard length or anal fin length as a covariate. The detected QTL were similar to those detected without covariates (Table S2 and Table S3), suggesting that these QTL were controlling the numbers of papillary processes even after correction for the standard length or the anal fin length.

Tests of QTL repeatability in the AFOM family
To test the QTL repeatability, we conducted QTL mapping on an independent F 2 family: the AFOM family. We found significant QTL on LG4, LG7, LG12, LG22, and LG23 ( Figure S3), which are summarized in Table 2. All QTL except those on LG4 and LG23 had the expected directions of effects (i.e., the southern alleles increased the number of papillary processes) ( Figure 6).
Only one QTL on LG22 overlapped with the QTL found in the OFAM family. A comparison of the heat maps of the LOD scores between the OFAM and AFOM families showed that few QTL were shared between families ( Figure 5). When we reduced the LOD significance threshold from P = 0.05 to P = 0.1, several suggestive QTL were found (Table S4): the QTL on LG22 (40 cM) for Ray 12, the QTL on LG22 (40 cM) for Ray 13, the QTL on LG7 (45 cM) for Ray 14, and the QTL on LG19 (20 cM) for Ray 15. Therefore, in addition to the QTL on LG22, the QTL on LG19 may also be shared between the two families.
Next, we tested whether the shared QTL on LG19 and LG22 had similar effects between families. The QTL on LG22 had a similar effect: in both families, it had no dominance effects (|d/a| , 0.5), and the southern allele increased the number of papillary processes (Table 1 and  Table 2). In contrast, the QTL on LG19 had different directions of effects: it had the opposite direction of effect in AFOM (Table S4), whereas it had the expected direction of effect in OFAM (Table 1).
In the AFOM family, we also found significant correlations between the papillary process number and the standard length or the anal fin length (Table S5). However, the inclusion of the standard length or the anal fin length as a covariate did not qualitatively change our results (Table S6 and Table S7).
Little overlap between QTL for papillary process number and anal fin length As described previously, there are phenotypic correlations between the raw values of the papillary process number and the anal fin length ( Figure S4). Because both the papillary process number and the anal fin length were correlated with the standard length, we calculated the residuals of these traits against the standard length. Significant correlations were found even between the residuals after we corrected for the standard length ( Figure S5), suggesting that some QTL for the papillary process number and the anal fin length might be shared.
In OFAM, we compared QTL for variations in papillary processes with previously identified QTL for anal fin length . The QTL on LG14 and LG19 overlapped (Figure 3), whereas other QTL on LG11, LG17, and LG22 did not. Overlapping QTL on both LG14 and LG19 had the expected directions of effects on both traits: the southern alleles increased both the anal fin length ) and the number of papillary processes (Figure 4).
Because we did not perform QTL mapping of the anal fin length and anal fin growth curves in the AFOM family in our previous study , we conducted QTL mapping of these traits in the AFOM. When the anal fin length was analyzed separately at each age, we found significant QTL on LG10 and LG15 (Table S8; Figure S6). When the weights of the orthogonal polynomial growth curves were analyzed, we found significant QTL on LG13 and LG14 in addition to LG10 and LG15 (Table S8; Figure S6). The QTL on LG14 and LG15 had the expected directions of effects (i.e., the southern alleles increased the anal fin length), whereas the QTL on LG10 and LG13 had the opposite directions of effects (i.e., the northern alleles increased the fin length) ( Figure S7). A comparison with QTL found in the OFAM family revealed that the QTL on LG10 and LG14 were shared between the AFOM and OFAM families. The directions of effects of the QTL on LG10 and LG14 were the same between the AFOM and OFAM families. Finally, we found no overlapping QTL between anal fin length/ growth curves and papillary processes in the AFOM family.

DISCUSSION
Independent genetic control of anal fin length and the number of papillary processes A major finding in our study was that there was little QTL overlap between anal fin length and the number of papillary processes in Japanese medaka. Thus, the two courtship traits are not necessarily controlled by the same genetic loci. Furthermore, we found that the numbers of papillary processes on different fin rays were controlled by different QTL. These results suggest that the different courtship traits associated with the anal fins can evolve independently. Consistent with this hypothesis, the genus Oryzias exhibits diversity in anal fins. Indonesian Oryzias species exhibit sexual dimorphism in anal fin length, with males having longer fins than females but lacking papillary processes (Iwamatsu 2006). Furthermore, Oryzias mekongensis males have less sexually dimorphic and shorter anal fins, although they still have papillary processes on the anal fins (Iwamatsu n  2006). The independent genetic control of each trait may enable the genus Oryzias to exhibit sexual dimorphism diversity in anal fin morphology.
The clustering of genes controlling functionally related traits may be favored when two divergent populations are faced with gene flow because physical linkage can prevent functionally related genes from being separated by recombination Jiggins and Bridle 2004;Yeaman and Whitlock 2011). For example, clusters of genes involved in reproductive isolation are found in linkage maps of sympatric pairs but not allopatric pairs Brown et al. 2004). In the present study, we used allopatric populations near the northern and southern limits of medaka distribution. Unless there was gene flow between the divergent populations, the physical linkage of QTL for different components of courtship would not necessarily be favored (Feder et al. 2013).

Repeatability of QTL
We noted that some QTL were found in two independent crosses whereas others were not. There are several possible reasons why QTL often fail to replicate across different families (Mauricio 2001;Barnwell and Noor 2008). First, grandparents are simply random samples from the parental populations and do not reflect all of the genetic variation existing in the parental populations (Barnwell and Noor 2008). Unlike inbred laboratory animals, many alleles with different effects are likely segregating within wild populations (Mauricio 2001). This may be particularly true when the analyzed traits are polygenic and/or selection is weak. Second, small sample sizes often lead to the failure of detection of significant QTL (Beavis 1994(Beavis , 1998. Although our sample sizes were relatively small, this factor alone could not explain the failure to replicate QTL because Figure 6 Comparison of the number of papillary processes among genotypes at the single-nucleotide polymorphism marker nearest the QTL in the AFOM family: SS, homozygotes of the southern alleles; SN, heterozygotes; NN, homozygotes of the northern alleles. n we did not find any suggestive QTL in one cross that overlapped significant QTL in the other cross, even when the significance threshold of the LOD score was reduced. Moreover, a comparison of the heat maps showed that the LOD score profiles apparently differed between the two families. Third, some epistatic interaction or environmental effects may mask the QTL effects (Vikram et al. 2011). OFAM is derived from a grandmother from Okinawa and a grandfather from Aomori, while AFOM is derived from a grandmother from Aomori and a grandfather from Okinawa. Therefore, if sex chromosomes (i.e., Y chromosomes) or maternal factors (i.e., mitochondrial DNA) had any epistatic interaction with QTL, we would detect different QTL between the OFAM and AFOM families. Consistent with this idea, Figure S4 indicates that the OFAM males (upper panel) tended to have a relatively smaller number of papillary processes than the AFOM males (lower panel). Although some QTL were not repeatable, we found several replicated QTL between families. Divergent alleles may segregate at relatively higher frequencies at these loci. It is also possible that the effects of these QTL may be independent of any epistatic or environmental effects.

Candidate genes
Several genes are known to control fin development or regeneration within 95% Bayesian credible intervals of QTL (Table 3). One QTL for the papillary process number that was replicated in both families contained bmp4. Inhibition of bmp4 has been suggested to disrupt fin growth during zebrafish fin regeneration (Smith et al. 2006). One QTL associated with both anal fin length and papillary process number in the OFAM family contained frem2a (Carney et al. 2010) and ar (Ogino et al. 2013). The ar gene encodes the androgen receptor (Ogino et al. 2013) and androgens can induce the expression of papillary processes on the anal fins of Japanese medaka (Ogino et al. 2013;Nakamura et al. 2014). Finally, a previous report demonstrated that variation in anal fin shape between a Japanese population and a Korean population of Oryzias latipes was associated with the genotype at cyp1b1 (Katsumura et al. 2014), and we found a significant QTL for anal fin length in the AFOM family near cyp1b1 on LG15 (Table S8). These genes would be good starting points for further investigation into the molecular mechanisms underlying the diversification of anal fin processes and length.
However, it should be noted that we only used 384 SNP panels in this study; therefore, we may have overlooked several QTL with small effects. Even if all of the QTL were captured, their precise locations could not be accurately mapped. Increasing the number of markers using next-generation sequencing technologies would help to identify the causative genes in future studies (Baird et al. 2008;Elshire et al. 2011;Peterson et al. 2012;Glazer et al. 2015;Kirchmaier et al. 2015). Elongation of fins can be achieved by increasing the number of segments and/or elongating each segment (Iovine 2007;Perathoner et al. 2014). Elongation of anal fins in the southern population of the Japanese medaka is mainly caused by an increase in the number of segments (Kawajiri and Yamahira 2011); therefore, further molecular developmental studies would also help to identify the causative genes. Due to the presence of phenotypic diversity between closely related species (Iwamatsu 2006;Mokodongan and Yamahira 2015) and the availability of gene manipulation technologies (Kinoshita et al. 2009;Ansai et al. 2013;Ansai and Kinoshita 2014), the genus Oryzias will serve as a good model system for further investigation into the molecular mechanisms underlying natural variation in sexually dimorphic or sex-specific traits between closely related species and populations.