1 INTRODUCTION

As an important commercial marine species, Japanese flounder (Paralichthys olivaceus) is widely cultured in China, Japan and Korea. However, due to accumulating stress from biological, physical and chemical factors during the aquaculture process, it regularly suffers from a disease that inflicts heavy, annual financial losses on producers. Edwardsiellosis, one of the most deadly diseases of Japanese flounder, is caused by a Gram-negative bacterium, Edwardsiella tarda, which was first recorded from Japanese eel by Hoshina (1962). Infected fish eventually develop ascitesosis (Zhang et al., 2003), leading to death. Despite its devastating effect on aquaculture, possible defense mechanisms and variations in disease resistance are, as yet, unknown, which has restricted the breeding of a disease-resistant strain.

Intensive efforts have been devoted to the rapid identification of markers linked with specific regions of the genome. Bulked segregant analysis (BSA) is emerging as a powerful tool to fulfill this role. The method is based on populations showing segregation in particular alleles, separating them into different gene pools according to phenotypes. Compared with traditional association analysis, its efficiency and simplicity has made it very popular in QTL mapping and association analysis (Mackay and Caligari, 2000).

In recent years, BSA-QTL has been widely used to detect various economic traits in farmed fish species, with significant progress being made on discovering markers related to phenotypic characteristics such as sex, heat tolerance, disease resistance and growth. In the northern snakehead (Channa argus), Liu et al.(2011) found a sex-related simple sequence repeat (SSR) marker that only appeared in females and hoped to culture pure female hybrids using femalespecific DNA fragments. Lu et al.(2007), combining BSA, screened nine heat-related molecular markers in P. olivaceus, three loci of which were significantly negatively correlated with upper thermal tolerance and six loci were positively correlated. Another heatrelated molecular marker was detected in turbot, Scophthalmus maximus (Ma et al., 2011), which may facilitate the selective breeding of high temperatureresistant strains. In the yellow croaker (Larimichthys crocea), Liu et al.(2013) found two molecular markers closely related to growth traits, and showed that both markers were located in LG11, which is a candidate linkage group with important genes related to growth traits. Xu et al.(2009) uncovered a microsatellite marker positively related to growth in S. maximus, which could use to select stocks for growth performance at a molecular level. With reference to disease resistance, Fuji et al.(2006) and Fan et al.(2014) identified several markers for resistance to lymphocystis virus disease (LD) in Japanese flounder. Fuji et al.(2006) considered a marker on LG15 as a candidate locus for markerassisted selective breeding to improve LD resistance. Wang et al.(2014) observed two markers located on LG7 associated with resistance to the fish pathogen, Vibrio anguillarum, in P. olivaceus ; a QTL candidate region was also identified but specific genes were not scanned. However, to date, there have been no studies on SSR markers or QTL related to E. tarda resistance in Japanese flounder, which is hampering the development of a disease-resistant strain.

The present study was designed by combining QTL mapping with BSA to screen disease-resistant markers and loci in a family challenged with E. tarda in 2013. The selection of families of Japanese flounder resistant to V. anguillarum was initiated in 2007 and 2009, but families specifically resistant to E. tarda were not selected until 2013(Zhang et al., 2003) and were based on families resistant to V. anguillarum (Wang et al., 2014). In that year, 32 of 56 Japanese flounder families were infected with E. tarda by intraperitoneal injection and those with different resistant capability revealed (Fig. 1). It was hoped that the application of these markers would improve Japanese flounder culture and facilitate the selection of an edwardsiellosis-resistant strain.

2 MATERIAL AND METHOD 2.1 Family selection and phenotype analysis

Four populations of Japanese flounder, e.g. Japanese stock, Korean stock, Chinese stock (resistant of V. anguillarum) and Yellow Sea wild stock, were used to established different families. Families established in 2007 were challenged with V. anguillarum and, of these families, F0768 had the highest survival rate (60%). In 2009, families with maximum survival abilities were considered a fifth population, i.e., after V. anguillarum infection, F0915 was selected as having the highest survival rate (72.82%). F1333, which was selected for screening markers for resistance to both E. tarda and V. anguillarum, was a cross of F0768(♀) and F0915(♂), and was expected. A total of 106 juveniles from F1333 were challenged with E. tarda and a further 31 families were established and challenged in the same year as controls.

Six-month-old juveniles, ~10 cm in length, were cultured in a 0.28-m³ tank under flow-through conditions with a constant supply of fresh water at 19±1℃. E. tarda was injected intraperitoneally at a concentration of 3.69×10⁵ CFU/mL, and the number and total length of dead fry recorded every 2 h after infection until no further mortalities were observed.

Two models were adopted for QTL mapping (Table 1). In model 1, the 47 dead juveniles were designated as 1 and survivors as 0(Houston et al., 2008 ; Ozaki et al., 2010). In model 2, the 10 dead individuals recorded within 48 h were designated as 4, the 25 dead juveniles observed at 48-96 h were designated as 3, the eight dead juveniles observed at 96-144 h were designated 2, the four dead individuals observed at 144-192 h were classified as 1, and all survivors were classified as 0(Shoba et al., 2012 ; Wang et al., 2014).

2.2 SSR marker selection and establishment of BSA

According to the genetic linkage map of Song et al.(2012), 180 SSR markers with an average interval of 10 cM were used in this study; the corresponding microsatellite primers were synthesized by Beijing Genomics Institute. Genomic DNA from flounders was extracted by the classical phenol-chloroform method (Liu et al., 2005). Resistant bulks (RB) were generated by combining equal volumes of DNA from 15 survivors at a concentration of 50 ng/mL. Likewise, susceptible bulks (SB) were generated from the 15 individuals that died first.

2.3 PCR amplification and PAGE

Polymerase chain reaction (PCR) amplification was performed in a 96-well plate at a 16-μL reaction volume containing 0.8 μL dNTP, 1.6 μL 10× buffer, 11.65 μL dd H₂ O, 0.15 μL Taq, 0.4 μL of each primer and 1 μL DNA at a concentration of 50 ng/μL. The program was 95℃ for 10 min for initial denaturation, followed by 30 amplification cycles at 95℃ for 30 s, 59.5℃ for 30 s and 72℃ for 30 s, then 5 min at 72℃ for the final extension. PCR products were stored at 4℃ for further use. For polyacrylamide gel electrophoresis (PAGE), 4.0 μL of each denatured sample was loaded on 6% denaturing PAGE gels and electrophoresed under a 1 500-V constant voltage for 1.5 h (Wang et al., 2014). After electrophoresis, silver staining was performed and the band patterns analyzed.

2.4 Verification

When a band appeared in a specific pool, it was regard as a disease-related band pattern. The markers yielding disease-related band patterns were re-tested on all individuals in the same pool and those failing to reproduce the pattern were excluded: positive markers were confirmed in 106 individuals. Fisher's exact chisquare tests for the different alleles between RB and SB were performed with cross tabulation using SPSS 17.0 software (chi-square test 95% confidence interval of 0.05)(SPSS Inc., Chicago, IL).

2.5 QTL mapping and data analysis

A genome scan was performed using the abovementioned 180 SSR markers in 30 animals in RB and SB, as well as their parents. Markers that displayed significant differences (via chi-square test) were genotyped in 106 challenged juveniles from F1333 and their parents - a total of 108 individuals. QTL analyses were performed using the two models described above. The logarithm of odds (LOD) threshold was fixed at 3.0, and the significant between marker and trait was set at P<0.05. Genotyped markers were mapped by JoinMap 4.0(Van Ooijen, 2006), and WinQTLCart 2.5 software (Wang et al., 2012) was used to analyze the QTLs associated with disease resistance.

3 RESULT 3.1 SSR and BSA analysis

A total of 180 SSR markers (Supplementary Table) covering 24 linkage groups in Japanese flounder (Song et al., 2012) were used to identify polymorphic bands between the RBs and SBs. After scanning, six SSR markers, i.e., Scaffold 270-13812, Scaffold 404- 21594 and Scaffold v404-21589 in LG1, Scaffold 410-1885 in LG15, Scaffold 826-5003 and Scaffold 364-18000 in LG8, were detected (Fig. 2). This result suggests that these six markers may be associated with disease resistance.

3.2 Verification

To verify differences in microsatellite loci between the RBs and SBs, 30 individuals from the two bulks were genotyped (Fig. 3). The SSR alleles were named A, B, C and D, according to their decreasing length. For example, the longest allele of Scaffold 404-21594 was named Scaffold 404-21594-A, while the shortest was designated Scaffold 410-1885-D. Different alleles were tested via chi-square tests in SPSS 17.0; SSR markers showing differences between the two bulks are listed in Table 2. From the chi-square test results, Scaffold 410-1885-D (P =0.067 9), Scaffold 826-5003(P =0.067 9) and Scaffold 364-18000(P =0.067 9) displayed no clear differences (P>0.05) between the two bulks, while Scaffold 404-21594(P =0.028 1), Scaffold 404-21589(P =0.009 9) and Scaffold 270-13812(P =0.028 1) showed significant differences (P<0.05). The frequency of Scaffold 404- 21589-B in the RB was 80% compared to 33.3% in the SB. Scaffold 404-21594-C and Scaffold 270- 13812-B had the same frequency (RB at 73.3% versus SB at 33.3%): their frequency in RB was significantly higher than in SB (P =0.028 1, <0.05) which means the alleles were more associated with resistant character than susceptible one. In contrast, Scaffold 404-21594-D and Scaffold 270-13812-A displayed the opposite pattern (RB at 26.7% versus SB at 66.7%), where the frequency in RB was significantly lower (P =0.028 1, <0.05). As further identification was required for all individuals, 106 challenged juveniles were genotyped. With the exception of Scaffold 410-1885-B, the other three markers were still significantly different (Fig. 4). For example, of 59 resistant juveniles, Scaffold 404-21594-C was found to have different bands in 42 individuals at a frequency of 71.2%; whereas 15 of 47 susceptible juveniles had different bands at a frequency of 32%. All three markers showed significant differences between the two bulks (P<0.001)(Table 3). Therefore, it was concluded that these three markers in LG1 may be closely associated with resistance to E. tarda infection but had no relationship with V. anguillarum in LG7(Wang et al., 2014)

3.3 Single marker and QTL map analysis

All 43 SSR markers in LG1 were genotyped in the 106 challenged juveniles because the three candidate markers (Scaffold 404-21594, Scaffold 404-21589 and Scaffold 270-13812) were located in LG1, and 24 of them were used to construct a new, small linkage map (Fig. 5). The length of the new linkage map was 92.4 cM. There was only one significant marker at the 5% and 0.01% levels, while there were seven and nine significant markers at the 1% and 0.1% levels, respectively, in model 1. These 18 markers were responsible for 0-15.75% of the phenotypic variance (R²) and the additive effect ranged from -0.379 26 to 0.413 398. In model 2, there were eight and nine markers at the 5% and 1% levels, respectively, accounting for 0.01%-10.35% of the phenotypic variance (R²) and from -1.224 71 to 1.252 39 of the additive effect. By comparing the two models, the highest R² value and additive effect were found at Scaffold 790-4911(Table 4).

Composite interval mapping (CIM) was adopted to re-analyze the association between the genotyped markers and disease resistance. A QTL with a LOD ≥3 was considered a promising candidate. Altogether, six QTLs were identified on LG1. Of these six QTLs, qE-1, qE-2 and qE-3 constructed in model 1, explained 74.3%, 17.0% and 77.1% of R², respectively, with additive effects of -0.997 43, -0.414 84 and 0.997 465, respectively, and peaked at 36.1, 66.7 and 86.5 cM with LOD scores of 58.4, 3.5 and 55.7, respectively. The QTLs - qE-4, qE-5 and qE-6 - constructed in model 2 explained 16.0%, 89.5% and 77.4% of R², respectively, with additive effects of 1.252 4, -2.968 and 2.894, respectively, and peaked at 66.7, 71.3 and 77.3 cM with LOD scores of 3.4, 21.5 and 13.7, respectively (Table 5 and Fig. 6). The phenotypic variance (R²) of the SSR markers located adjacent to the QTLs ranged from 0.07% to 16.91%. Scaffold 404-21594 explained the highest level of R²(16.91%)(Table 6). Coincidentally, it was found that qE-2 in model 1 and qE-4 in model 2 were both located at the same position of 66.7 cM, close to Scaffold 790-4911. After gene scanning of the Japanese flounder genome, two immune genes, Gvin1(interferon-induced very large GTPase 1) inside qE-2 and qE-4, and IL-8(interleukin-8) inside qE-5 were identified.

4 DISCUSSION

A thorough understanding of the genetic basis of pathogenic disease mechanisms is essential to improve production of aquaculture fish species, such as Japanese flounder (Paralichthys olivaceus). The detection of candidate genes by QTL mapping could provide some clues on the pathogenicity of disease organisms. To reduce extraneous factors, in this study, 6-month-old juvenile fish fry were selected and cultured in the same pond. Bulks of phenotypic extremes were made by obtaining samples from surviving and dead individuals (Wang et al., 2014).

Currently, a large number of markers linked to quantitative traits have been found in fish species using BSA (Lee et al., 2004 ; Keyvanshokooh et al., 2007 ; Wang et al., 2009a, b ; Liu et al., 2011). The technique is cost-efficient, simple to use and timesaving, which will greatly advance research in this field. In our study, three markers - Scaffold 404- 21594, Scaffold 404-21589 and Scaffold 270-13812 in LG1 - potentially linked to a disease-resistance gene in Japanese flounder were identified using BSA combined to SSR and single-marker analysis. However, none of the markers was able of fully demonstrating disease-resistance, because disease traits are quantitative and always controlled by several rather than a single gene. The F1333 was established based on familial resistance to V. anguillarum, but the three screened markers were located on LG1, and not on LG7, where the markers associated with resistance of V. anguillarum are located (Wang et al., 2014). This indicated that these markers were only associated with resistance to E. tarda rather than E. tarda and V. anguillarum ; in other words, the identified markers had nothing to do with V. anguillarum. Despite some deviation, the overall trend indicated that individuals may have the resistance gene but are not phenotypically resistant to the disease, which is in agreement with the studies of Xu et al.(2009) and Fan et al.(2014). The most economically important quantitative traits are continuously changing because they were usually controlled by multiple genes and easily influenced by environmental factors, and there is no clear association between phenotype and genotype (Xu et al., 2009 ; Yue, 2014). Unfortunately, in this study, only one family displayed significant differences between RB and SB. Thus, further evidence is required to verify the accuracy and feasibility of the process in other populations (Geldermann, 1975 ; Wang et al., 2008 ; Moen et al., 2009). Once the presence of resistance gene markers is confirmed, they can be used in marker-assisted selective (MAS) breeding.

Two different models were used for phenotypic analysis in this study, as described by Wang et al.(2014). Model 1 is normally used to analyze quantitative traits controlled by a few major genes, while model 2 is always used to analyze traits controlled by multiple minor genes. Each model identified three QTLs. Both qE-2 and qE-4 were located at position 66.7 cM, suggesting that this area may be more highly linked to E. tarda disease resistance than other QTLs. Thus, this region should be the primary focus for screening disease-resistant candidate genes. However, the possibilities of immune genes being located in other QTLs should not ignore; nevertheless, this study failed to located genes in qE- 1, qE-3 or qE-6.

Gvin1, the prototype member of the newest family of very large IFN-inducible GTPases (VLIG)(Joyce et al., 2009), is a member of the GTPase family. It is regulated by IFN-γ-induced genes, and is involved in the survival of pathogens residing in phagosomes vacuoles (Singh et al., 2007). However, in Japanese flounder, the relationship between Gvin1 and E. tarda is unclear. In our study, this gene was located inside qE-2 and qE-4, which provides solid evidence to further uncover its novel role in E. tarda resistance. IL-8 is also a promising candidate gene because it is located in the qE-5 region in LG1 and was shown to be expressed in specific tissues, such as head kidney and spleen, after stimulation in Japanese flounder (Lee et al., 2001). Produced by peripheral blood mononuclear cells (Matsushima et al., 1988), macrophages (Goodman et al., 1991), lymphocytes (Gregory et al., 1988) and epithelial cells (Nakamura et al., 1991), the targets of IL-8 are also located in multiple cells. Similar to Gvin1, further detailed analysis of the relationship between IL-8 and E. tarda resistance in Japanese flounder is necessary.

In this study, the highest R² value was 15.75%, which fails to meet the requirements for MAS breeding. In addition, the number of microsatellite markers is insufficient and improved high-density, integrated genetic maps are urgently required, which should include more families or populations (Yue, 2014), more accurate phenotype collation data and better positioning methods. As emerging markers, single nucleotide polymorphisms (SNPs) are more widely distributed and are more tightly linked to genes compared with SSRs. Using next-generation sequencing (NGS) techniques, it is now possible to develop a large number of SNPs for linkage mapping (Miller et al., 2007, 2012 ; Davey et al., 2011). Genome-wide association studies (GWAS)(Hirschhorn et al., 2005) are currently being applied to analyze individual SNP variations within the whole genome. The process can detect minute effects of marker-trait associations, including the specific genes associated with disease; however, the high financial cost is limiting its application (Mccarthy et al., 2008). As an alternative, the method of combining QTL with GWAS research is relatively cost-effective and efficient, so may play a more significant in MAS breeding in the future.

5 CONCLUSION

In conclusion, three SSR markers associated with disease resistance to E. tarda in Japanese flounder were detected and reported for the first time. Future studies should focus on screening for disease-resistant genes located between the flanks of QTLs. Moreover, a high-density genetic linkage map (e.g. an SNP map) is required to identify more markers highly linked with E. tarda resistant QTLs. These development and advances would greatly facilitate the breeding of an edwardsiellosis-resistant cultivar of Japanese flounder.