Evolutionary Ecology Research, 2013, 15: 295–311

Are Japanese freshwater populations of
threespine stickleback derived fromthe Pacific Ocean lineage?
Lara M. Cassidy1,2, Mark Ravinet1,3, Seiichi Mori4 and Jun Kitano1,5
1

Ecological Genetics Laboratory, National Institute of Genetics, Mishima, Japan,
Trinity College, Dublin, Ireland, 3Queen’s University Belfast, Belfast, UK, 4Gifu-keizai University,
Ogaki, Japan and 5PRESTO, Japan Science and Technology Agency, Kawaguchi, Saitama, Japan

2

ABSTRACT
Background: The presence of ecological opportunity can trigger adaptive radiation. Freshwater colonization of marine ancestors during the post-glacial dispersal could have triggered
adaptive radiation in the threespine stickleback (Gasterosteus aculeatus). Japanese marine
threespine stickleback can be classified into two genetically divergent groups, the Pacific Oceangroup and the Japan Sea group.
Question: Are all Japanese freshwater threespine stickleback populations derived from the
Pacific Ocean lineage? Or are some freshwater populations derived from the evolutionarily
divergent Japan Sea lineage?
Methods: We collected stickleback from 22 different locations across Japan, including
nine freshwater populations, five Pacific Ocean anadromous populations, and eight Japan
Sea anadromous populations. We determined the genotypes for 11 different microsatellite
markers. We created phylogenies using different measures of genetic distance based on bothallele frequencies and allele lengths. We inferred population structure using Bayesian analysis.
Results: All freshwater populations analysed were genetically similar to the Pacific Ocean
anadromous populations, suggesting that they are likely derived from the Pacific Ocean lineage
rather than the Japan Sea lineage.
Keywords: adaptive radiation, microsatellites, phylogeography, stickleback.

INTRODUCTION
The presence of ecological opportunity can trigger adaptive radiation (Schluter, 2000; Glor, 2010;
Losos, 2010). When interspecific competition or predation is reduced, phenotypic diversifica-tion occurs as a result of colonization and adaptation to vacant niches. Newly formed lakes
and islands or mass extinction can provide ecological opportunities to many lineages(Schluter, 2000; Gillespie and Baldwin, 2010; Glor, 2010; Losos, 2010). However, not all lineages can exploit

Correspondence: J. Kitano, Ecological Genetics Laboratory, National Institute of Genetics, Yata 1111, Mishima,
Shizuoka 411–8540, Japan. e-mail: jkitano@lab.nig.ac.jp
Consult the copyright statement on the inside front cover for non-commercial copying policies.
© 2013 Jun Kitano

296

Cassidy et al.

these opportunities. For example, Darwin’s finches and Hawaiian honeycreepers haveachieved tremendous adaptive radiations, while Galápagos mockingbirds and Hawaiian
thrushes have not (Losos, 2010). Differences in dispersal rates or geographical proximity tohabitats with vacant niches play an important role in shaping colonization history. Lineages
that are able to colonize novel environments and diversify rapidly may therefore preclude
later colonization by other competing lineages (Glor, 2010; Losos, 2010). Alternatively, colonizing
taxa might lack the key traits or genetic variation that allows exploitation of vacant niches(Glor, 2010; Losos, 2010). Interactions between adaptation and colonization success mean thatthe phylogeography of evolutionary divergent lineages is of considerable importance(Waters, 2011). Comparisons between lineages that successfully diversified and those that did
not are essential for a better understanding of the factors promoting and constraining
adaptive radiations.
Temperature shifts, ice sheet movement, and sea-level change during the Quaternary
period altered the range and distribution of taxa, thus influencing gene flow, genetic
diversity, and divergence (Avise, 2000; Bennett, 2008). The glacial and inter-glacial cycles of the
Pleistocene have driven adaptive radiations in a number of post-glacial temperate northernfish species (Schluter, 1996). Ice retreat following the Last Glacial Maximum exposed numerous
freshwater habitats that were colonized by multiple ancestral marine fish species, providing
novel niches and leading to rapid evolutionary diversification (Taylor, 1999; Schluter, 2000). One of
the best studied is the threespine stickleback (Gasterosteus aculeatus L.) species complex,
characterized by independent colonization of freshwater environments by marine ancestralforms throughout its circumpolar distribution, leading to the parallel evolution of extensive
phenotypic diversity (Wootton, 1976; Bell and Foster, 1994; McKinnon and Rundle, 2002). Althoughmarine ancestral sticklebacks are relatively homogeneous in phenotypic traits, freshwater
sticklebacks diversified and adapted to diverse environments after colonization of a varietyof newly created freshwater environments (Wootton, 1976, 1984; Bell and Foster, 1994).
Japanese marine threespine stickleback can be classified into two genetically divergent
groups, the Pacific Ocean anadromous form (PA) and the Japan Sea anadromous form (JA)
(Higuchi and Goto, 1996; Kitano et al., 2007, 2009). Around the Japanese archipelago, the distribution of
JA fish centres primarily on the Sea of Japan and Sea of Okhotsk, with populations
found across the western seaboard of Honshu (Higuchi and Goto, 1996). In contrast, PA fish havea much broader distribution, occurring across the Pacific Ocean from North America toJapan (Higuchi and Goto, 1996). Both forms are sympatric in Hokkaido, the northern islands of
the Japanese archipelago, sometimes co-occurring in the same drainage system (Higuchi and
Goto, 1996; Kitano et al., 2007). PA and JA fish probably diverged during a period of extended
allopatry 0.5–2 million years ago when sea-level change isolated the Sea of Japan from the
Pacific Ocean (Yamada et al., 2001; Kitano et al., 2007).
In addition to the two anadromous forms, numerous freshwater populations of threespine
stickleback with divergent morphological and life-history traits also occur in the northernparts of the Japanese archipelago (Mori, 1987, 2003). However, the evolutionary history ofJapanese freshwater populations remains unclear. Previous attempts to resolve the phylogeography of Japanese stickleback populations have provided conflicting evidence of ancestry. Allozyme analysis demonstrated that two clearly distinct lineages occurred across theregion (Higuchi and Goto, 1996). A strong correlation between genetic data and morphologicalmeasurements (caudal plate height and caudal keel ossification) supports the hypothesisthat these two lineages correspond to the JA and PA forms. Furthermore, all of the freshwater populations surveyed in this study grouped clearly within the PA lineage (Higuchi andOrigin of Japanese freshwater stickleback populations

297

Goto, 1996). These data suggest that the PA rather than the JA lineage might be the ancestor of

the Japanese freshwater populations. However, a further study using mitochondrial DNA on
similar populations did not detect such clear genetic divergence (Yamada et al., 2001, 2007). Instead,all Japanese populations were grouped in a monophyletic lineage consistent with theJapanese clade revealed by a global survey of threespine stickleback mitochondrial DNA(Ortí et al., 1994). This suggests that introgression of mitochondrial haplotypes between PA
and JA has likely occurred on secondary contact (Yamada et al., 2001, 2007). Given conflicting
results from different marker types, further evidence from nuclear loci is required.
In light of the complicated evolutionary history of Japanese stickleback populations,
we investigated the genetic relationships between marine and freshwater populations
using microsatellite markers. In short, we wished to examine whether the Japanese freshwater populations are derived from the Pacific Ocean lineage, the divergent Japan Sea
lineage or both.
MATERIALS AND METHODS
Sampling and fish collection
Threespine stickleback were collected from 22 sampling locations across Japan (Fig. 1,
Table 1). This included nine freshwater habitats, eight JA stickleback habitats, and five
PA stickleback habitats. All freshwater populations analysed were sampled from habitatsthat were landlocked or located above dams and falls. JA and PA forms were collected fromthe estuary and first classified into JA and PA forms on the basis of external morphology:
PA fish have a larger caudal plate height and a larger body size than JA fish (Higuchi and Goto,
1996; Kitano et al., 2007). Ten individuals from each population were sampled, with the exceptionof the Gifu population from which only eight samples could be obtained. DNA was isolated
from the pectoral fin using the DNeasy Blood and Tissue Kit (Qiagen, Hilden, Germany) as
described previously (Adachi et al., 2012).
Microsatellite marker amplification
Fish were genotyped using 12 microsatellite markers (Stn170, Stn215, Stn233, Stn64, Stn76,
Stn159, Stn46, Stn90, Stn120, Stn278, Stn332, and Stn384), chosen because they arelocated on different threespine stickleback linkage groups (Peichel et al., 2001), are not linked
to sex (Peichel et al., 2004; Kitano et al., 2009), and were observed to be polymorphic in many
populations previously analysed (Kitano et al., 2008a, 2009; Adachi et al., 2012). Forward primers werelabelled with HEX, NED or FAM and the 5⬘-ends of the reverse primers were tailed with
GTTTCTT to increase the accuracy of fragment analysis (Ballard et al., 2002). Markers werecombined into three different multiplexes based on dye colour and were then amplified
using the KAPA2G Fast Multiplex PCR Kit (KAPA Biosystems, Woburn, MA, USA).
After 3 min of 95⬚C, 30 cycles of 95⬚C for 15 s, 60⬚C for 30 s, and 72⬚C for 30 s were
performed, followed by 10 min of 72⬚C. Amplified fragments were analysed by BEX Co.
Ltd. (Tokyo, Japan). Allele lengths were then determined and scored using Peak Scanner
Software (Life Technologies, Grand Island, NY, USA). Eleven of the 12 primer pairs
amplified microsatellite fragments in all populations, whereas one marker (Stn76) failed
to amplify in JA populations. Alleles at this locus were therefore private to PA and freshwater populations and thus were excluded from genetic diversity and distance analyses.

298

Cassidy et al.

Fig. 1. Map showing the sampling sites. Open circles and solid circles denote the allopatric JA and PA
sampling sites, respectively. Grey circles denote the sympatric PA/JA sampling sites. Squares denote
freshwater populations. In the upper panel, proposed distributions of the JA and PA lineages areshown in red and blue respectively. In the lower panel, the sampling points of Hokkaido Island andthe northern parts of Honshu Island are shown.

Origin of Japanese freshwater stickleback populations

299

Table 1. Information relating to the sampling sitessPopulation

Forms

LatitudeLongitude

Collection date

Sample size

AisakaAizu
Chimikeppu
Gensui
Gifu
Kussharo
Nishitappu
Shikotsu
TowadaHyotanWatarichippu
Choboshi
Akkeshi PA
Harutori PA
Abashiri
Benkei
Notojima
Ogawarako
Sarufutsu
Onnechikappu
Akkeshi JA
Harutori JA

Freshwater
Freshwater
Freshwater
Freshwater
Freshwater
Freshwater
Freshwater
Freshwater
Freshwater
Allopatric PA
Allopatric PA
Allopatric PA
Sympatric PA
Sympatric PA
Allopatric JA
Allopatric JA
Allopatric JA
Allopatric JA
Allopatric JA
Sympatric JA
Sympatric JA
Sympatric JA

40.59237.51143.62939.36535.39443.60042.64642.77540.44543.03243.036443.257943.10542.96943.96042.82537.12740.84045.25542.96643.06842.969141.221139.866
143.885
141.897136.620144.348141.476
141.400
140.842144.844145.0532
145.556144.891
144.396
144.200
140.188137.044141.372
142.239144.116144.884
144.396

October 2010
May 2008
June 2011
August 2011
May 2007
May 2003
September 2010
September 2010
June 2011
June 2010
June 2011
May 2012
June 2006
June 2006
May 2008
June 2008
May 2007
March 2002
June 2005
June 2001
June 2006
June 2010

10
10
10
10
8
10
10
10
10
10
10
10
10
10
10
10
10
10
10
10
10
10

However, the marker was used as prior information for STRUCTURE analysis of
population structure (see below).
Genetic diversity, distance, and population clustering
Populations were tested for linkage disequilibrium and Hardy-Weinberg equilibrium(HWE) using Genepop (http:genepop.curtin.edu.au/) (Raymond and Rousset, 1995). GenoDive
was used to estimate levels of genetic diversity (Ho , observed heterozygosity; He , expected
heterozygosity) and effective number of alleles (AE ) within populations (http:www.bentley
drummer.nl/software/Home.html) (Meirmans and Van Tienderen, 2004). To indicate genetic differentiation between populations, pairwise Weir and Cockerham’s (1984) FST and Jost’s D values(Jost, 2008) were estimated using the diveRsity R package (Keenan, 2012). General linear mixedmodels (GLMMs) were used to test for differences in AE , Ho , and He between freshwaterand marine populations. All statistical analyses were conducted in R 2.15.1 (R Development
Core Team, 2012).
To investigate population structure, we used a Bayesian clustering method implemented
in STRUCTURE (Pritchard et al., 2000). The number of assumed populations (K) was set from
2 to 23. The model allowed for admixture and prior information was provided, based on
private alleles at the Stn76 locus, which amplified in PA populations and freshwater
populations only. Without including the prior information, we obtained similar results. For
each run of the model, a burn-in of 50,000 steps was used, followed by 500,000 iterations.
For each value of K, the model was repeated 10 times to ensure consistent results. We

300

Cassidy et al.

evaluated the most probable value of K using log posterior probabilities and Evanno’s ∆Kmethod (Evanno et al., 2005) implemented in STRUCTURE Harvester (Earl and vonHoldt, 2012).
Independent runs were grouped using CLUMPP (Jakobsson and Rosenberg, 2007) and displayedusing DISTRUCT (Rosenberg, 2004).
To further examine the relationships between populations, we analysed the distributionof alleles across all loci using the program ADZE (Szpiech et al., 2008). Prior to analysis, populations were grouped as being Pacific Ocean, Japan Sea or freshwater. Following this, we used
ADZE to estimate the mean number of alleles per locus in each group and then the mean
number of private alleles. Unlike other rarefaction methods (Petit et al., 1998; Kalinowski, 2004), the
ADZE method calculates private alleles as those occurring private to all populations withina grouping. Furthermore, ADZE conducts this analysis at different sample sizes (g), therefore demonstrating the effect of sample size on these estimates; for the present analysis, thevalue of g was set to vary from 2 to 50. Finally, we used ADZE to estimate the mean numberof private alleles for all possible pairwise combinations of the Pacific Ocean, Japan Sea, andfreshwater groups. Following Szpiech et al. (2008), we hypothesized that more closely relatedgroups would retain a higher proportion of private alleles than those distantly related.
Microsatellite phylogenetic analysis
Two phylogenetic trees were generated, based on different measures of genetic distance. The
first was Nei’s DA (Nei et al., 1983), which is based solely on allele frequencies. The other was
stepwise weighted genetic distance (DSW) (Shriver et al., 1995), which takes into account theamplified microsatellite nucleotide length. DSW is considered more suitable for estimationof deeper phylogenetic split than Nei’s DA (Goldstein et al., 1995; Terazaki and Nei, 1996). Using
these two genetic distances, phylogenetic trees were created using neighbour-joining
methods (Saitou and Nei, 1987) with 3000 bootstrap replications in Poptree2 software(http:www.med.kagawa-u.ac.jp/∼genomelb/takezaki/poptree2/index.html) (Takezaki et al., 2010).
RESULTS
Genetic diversity and distance
All populations were in HWE for at least nine of the 11 loci following Bonferroni corrections (i.e. P > 0.005; see Appendix 1) and no linkage disequilibrium was detected between
any microsatellite loci (P > 0.05). Tests for genetic diversity showed that freshwater
stickleback populations have lower numbers of alleles and lower levels of heterozygosity
2
than marine populations (GLMMs with population as a random effect; AE : R = 0.73,
2
2
F1,19 = 35.81, P < 0.0001; Ho: R = 0.99, F1,19 = 33.82, P < 0.0001; He: R = 0.68, F1,19 =44.05, P < 0.0001; see Table 2). Within the marine populations, JA populations appeared
to have a higher mean number of alleles than their PA counterparts (mean AE ± ..:JA = 6.27 ± 0.83; PA = 3.63 ± 0.13), although this was not statistically significant (R 2 = 0.76,
F1,1 = 47.94, P = 0.09). In terms of heterozygosity, Ho did not differ between the two forms(P = 0.18) whereas He was higher in JA fish (R 2 = 0.83, F1,1 = 1524.5, P = 0.02).
Pairwise measures of genetic differentiation ranged from 0.00 to 0.70 (FST) and 0.95
(Jost’s D; see Table 3). The highest differentiation was present between freshwater
populations and between freshwater and marine populations. For example, differentiation
between JA fish from Akkeshi and the freshwater population in Aisaka was 0.30 and 0.82

Origin of Japanese freshwater stickleback populations

301

Table 2. Number of alleles (NA), effective number ofalleles (Effective NA), and expected heterozygosity (He)Population

NA

Effective NA

He

Freshwater
AisakaAizu
Chimikeppu
Gensui
Gifu
Kussharo
Nishitappu
Shikotsu
Towada
Mean

2.818
1.727
3.909
3.273
2.636
5.182
3.545
3.182
2.273
3.172

2.310
1.386
2.842
2.128
2.008
3.614
2.468
2.356
1.758
2.319

0.509
0.198
0.475
0.478
0.360
0.716
0.566
0.535
0.392
0.470

Japan Sea
Abashiri
Akkeshi JA
BenkeiHarutori JA
Notojima
Ogawarako
Onnechikappu
Sarufutsu
Mean

9.000
8.700
8.400
10.00
7.000
7.600
9.600
8.700
8.625

6.344
6.352
6.220
7.195
4.598
5.830
7.259
6.377
6.272

0.849
0.841
0.852
0.860
0.804
0.858
0.872
0.847
0.848

Pacific Ocean
Akkeshi PA
ChoboshiHarutori PA
HyotanWatarichippu
Mean

6.182
6.000
5.909
5.727
5.000
5.764

3.642
3.726
3.590
3.774
3.421
3.631

0.713
0.677
0.725
0.722
0.715
0.710

(FST and Jost’s D, respectively). In contrast, much lower differentiation occurred between JA
populations (Table 3). Between groups of populations, genetic differentiation was highest
between JA and freshwater populations (FST = 0.18, Jost’s D = 0.66) and between pooled
JA and PA populations (FST = 0.17, Jost’s D = 0.65). For both measures, divergence wasmarkedly lower between freshwater and PA populations (FST = 0.04, Jost’s D = 0.08).
Population structure
K = 4 was determined as the most probable number of populations by using the mean log
P(K), Evanno’s ∆K, and the q-values at each probable K (see Fig. 2). ∆K values indicatedthat K = 10 might be the true K, and several freshwater populations started to appear as
independent clusters. At both K = 4 and K = 10, JA and PA fish formed two distinct clusters.
Because we wanted to determine whether the freshwater populations are geneticallysimilar to the JA or PA, we analysed the genetic structure at K = 2. At K = 2, two clusters

0.769
0.809
0.759
0.732
0.828
0.633
0.741
0.796
0.777
0.002
0.000
0.605
0.415
0.000
0.670
0.605
0.000
0.659
0.000
0.144
0.015
—

SHI TOW ABA AKJ AKP CHO HRJ HRP HYO SAR WAT BEN NOT OGA ONN— 0.484 0.402 0.175 0.811 0.186 0.239 0.201 0.034 0.918 0.824 0.166 0.196 0.757 0.116 0.129 0.805 0.191 0.857 0.709 0.891
0.504 — 0.573 0.364 0.788 0.429 0.491 0.541 0.482 0.902 0.852 0.532 0.459 0.924 0.524 0.520 0.910 0.592 0.911 0.946 0.937
0.371 0.563 — 0.398 0.924 0.264 0.340 0.255 0.560 0.695 0.781 0.304 0.306 0.726 0.234 0.384 0.715 0.319 0.886 0.834 0.832
0.219 0.475 0.388 — 0.755 0.081 0.152 0.165 0.213 0.821 0.779 0.168 0.165 0.812 0.210 0.240 0.788 0.180 0.785 0.603 0.823
0.513 0.671 0.566 0.520 — 0.753 0.770 0.767 0.842 0.826 0.843 0.759 0.775 0.833 0.832 0.813 0.724 0.850 0.765 0.896 0.801
0.140 0.358 0.243 0.102 0.370 — 0.079 0.146 0.277 0.776 0.682 0.082 0.023 0.730 0.071 0.064 0.734 0.048 0.762 0.708 0.807
0.211 0.463 0.272 0.192 0.470 0.065 — 0.052 0.329 0.796 0.766 0.121 0.063 0.803 0.274 0.216 0.781 0.171 0.861 0.711 0.874
0.200 0.482 0.250 0.196 0.481 0.102 0.070 — 0.315 0.807 0.831 0.187 0.176 0.719 0.236 0.238 0.735 0.243 0.832 0.590 0.822
0.082 0.609 0.456 0.311 0.598 0.228 0.330 0.330 — 0.934 0.774 0.238 0.281 0.775 0.225 0.226 0.802 0.280 0.889 0.843 0.944
0.303 0.440 0.312 0.307 0.344 0.179 0.253 0.266 0.373 — 0.001 0.661 0.687 0.000 0.757 0.756 0.016 0.730 0.001 0.206 0.000
0.295 0.442 0.316 0.303 0.346 0.176 0.257 0.266 0.352 0.000 — 0.557 0.573 0.000 0.726 0.669 0.005 0.608 0.003 0.134 0.001
0.127 0.389 0.216 0.137 0.376 0.041 0.095 0.117 0.210 0.184 0.178 — 0.004 0.642 0.009 0.002 0.705 0.001 0.751 0.733 0.742
0.157 0.390 0.227 0.135 0.407 0.021 0.060 0.117 0.247 0.199 0.199 0.020 — 0.632 0.049 0.014 0.563 0.085 0.723 0.667 0.655
0.282 0.437 0.305 0.296 0.339 0.169 0.250 0.250 0.343 0.000 0.000 0.171 0.193 — 0.723 0.618 0.001 0.709 0.002 0.149 0.016
0.132 0.386 0.192 0.136 0.388 0.057 0.138 0.138 0.223 0.185 0.182 0.012 0.049 0.176 — 0.020 0.785 0.011 0.737 0.607 0.794
0.105 0.389 0.248 0.162 0.390 0.048 0.116 0.136 0.215 0.187 0.186 0.014 0.041 0.174 0.041 — 0.736 0.051 0.653 0.630 0.741
0.284 0.443 0.309 0.300 0.334 0.171 0.248 0.251 0.344 0.008 0.002 0.178 0.190 0.000 0.183 0.181 — 0.699 0.043 0.135 0.005
0.164 0.416 0.220 0.140 0.399 0.044 0.130 0.130 0.231 0.187 0.177 0.013 0.056 0.175 0.019 0.051 0.175 — 0.781 0.677 0.807
0.294 0.454 0.329 0.309 0.337 0.174 0.262 0.266 0.367 0.001 0.000 0.184 0.203 0.008 0.183 0.180 0.005 0.184 — 0.184 0.000
0.296 0.471 0.337 0.297 0.376 0.187 0.265 0.249 0.376 0.048 0.045 0.197 0.211 0.047 0.188 0.195 0.044 0.188 0.047 — 0.123
0.304 0.463 0.333 0.313 0.348 0.184 0.268 0.271 0.380 0.006 0.009 0.189 0.206 0.024 0.185 0.189 0.015 0.190 0.000 0.036 —
0.264 0.421 0.295 0.273 0.325 0.147 0.226 0.239 0.330 0.000 0.000 0.150 0.162 0.000 0.153 0.152 0.000 0.156 0.000 0.035 0.007

AIZ CHI GEN GIF KUS NIS

Note: Population abbreviations are as follows: AIS = Aisaka, AIZ = Aizu, CHI = Chimikeppu, GEN = Gensui, GIF = Gifu, KUS = Kussharo, NIS = Nishitappu, SHI =Shikotsu, TOW = Towada, ABA = Abashiri, AKJ = Akkeshi JA, AKP = Akkeshi PA, CHO = Choboshi, HRJ = Harutori JA, HRP = Harutori PA, HYO = Hyotan,
SAR = Sarufutsu, WAT = Watarichippu, BEN = Benkei, NOT = Notojima, OGA = Ogawarako, ONN = Onnechikappu.

AIS
AIZ
CHI
GEN
GIF
KUS
NIS
SHI
TOW
ABA
AKJ
AKP
CHO
HRJ
HRP
HYO
SAR
WAT
BEN
NOT
OGA
ONN

AIS

Table 3. Pairwise FST and Jost’s D values between populations (values below the diagonal are Weir and Cockerham’s (1984) estimated FST; values above thediagonal are estimated Jost’s D)

Fig. 2. Results of STRUCTURE analysis. Analysis of posterior probabilities with mean log P(K) ± .. (A) and ∆K (B) indicating that K = 4 orK = 10 is the true K; STRUCTURE output for K = 2, i.e. JA and PA lineages (C), K = 4 (D), and K = 10 (E).

304

Cassidy et al.

represented grouping into the JA and PA populations (Fig. 2C). Interestingly, all freshwater
populations, with the exception of the Gifu population, grouped with the PA cluster.
In contrast, all individuals from Gifu clustered with JA fish.
Shared private alleles between forms
Analysis of allelic distribution revealed that the number of mean distinct alleles per locus( ± ..) was greater for JA populations (14.53 ± 2.05) than the pooled PA or freshwater
populations (8.90 ± 1.41 and 8.22 ± 1.32 respectively; Fig. 3A). A similar pattern was seenfor private alleles (Fig. 3B). As expected, both measures of allelic distribution increased as afunction of sample size, although this was less pronounced for private alleles shared
between populations within a group (Fig. 3B). Tests of all pairwise combinations of
populations revealed the highest proportion of shared private alleles was between PA and
freshwater populations, while the lowest proportion occurred between JA and freshwater
populations (see Table 4).

Fig. 3. Results of Allelic Diversity AnalyZEr (ADZE) analysis for (A) mean distinct alleles per locus
and (B) mean private alleles per locus, each with increasing sample size (g). Dark grey circles denotepooled Japan Sea populations, light grey circles indicate Pacific Ocean populations, and open circles
indicate freshwater populations.
Table 4. Mean number of private alleles per locus for all pairwise combinationsof Japan Sea, Pacific Ocean, and freshwater populations (g = 50)
Combination
Pacific Ocean (PA) and freshwater
Pacific Ocean (PA) and Japan Sea (JA)Japan Sea (JA) and freshwater

Mean private alleles per locus

..

2.72
1.10
0.86

2.38
1.89
0.83

Origin of Japanese freshwater stickleback populations

305

Microsatellite phylogeny
In the tree of Nei’s DA (Fig. 4A), JA populations were monophyletic. All of the freshwater
populations were genetically similar to the PA populations. The tree also indicates that a
freshwater population from Gifu was located close to the root of the Japanese threespine
stickleback populations, suggesting that this Gifu population might have diverged earlier
than the other freshwater populations. In the tree of DSW (Fig. 4B), JA populations wereagain monophyletic and all freshwater populations were clustered with the PA populations.
The DSW tree also demonstrated that the Gifu population might have diverged earlier thanthe other freshwater populations. Furthermore, both trees do not indicate a monophyletic
origin of freshwater populations, but rather support multiple colonizations of fresh water in
the Japanese Archipelago.

DISCUSSION
Our microsatellite data confirm that the Japan Sea and Pacific Ocean anadromous stickleback forms are genetically divergent across their distribution, as previously suggested byallozyme analysis (Higuchi and Goto, 1996). Therefore, the discrepancy between nuclear DNA
data and mitochondrial DNA data is likely due to the introgression of mitochondrial DNA.
Despite the fact that this suggests gene flow throughout the evolutionary history of thisspecies pair, reproductive isolation persists where the forms occur in sympatry. At least five
types of sympatric/parapatric species pairs at varying stages of divergence have evolved
throughout the stickleback distribution (Hendry et al., 2009). Of these, the JA–PA species pairlies at the far end of the stickleback speciation continuum. The Japanese species pair is
unique in that it is the only known example of near-complete reproductive isolation anddivergence between anadromous populations within the stickleback species complex(McKinnon and Rundle, 2002; Hendry et al., 2009). Hybrid males produced in crosses between JA females
and PA males are sterile with attenuated sperm production (Kitano et al., 2007). The two forms
also exhibit divergent courtship behaviours, which act as a prezygotic isolating barrier when
the forms occur in sympatry (Kitano et al., 2007, 2008b, 2009). In addition, JA and PA sticklebacks
are divergent in body size and ecologically relevant traits such as number of gill rakers, diet,spawning habitat choice, and migration patterns (Kume et al., 2005, 2010; Kitano et al., 2007, 2009).
Divergent ecological selection between different spawning habitats and foraging regimes
appears to exist between them (Kume et al., 2010; M. Ravinet and J. Kitano, unpublished data). Therefore,
a combination of these multiple isolating barriers likely prevents these species fromhybridization in sympatry (Kitano et al., 2009).
We also found that all freshwater populations analysed were genetically similar to the PA
populations, suggesting that they are likely derived from the PA lineage rather than theJA lineage. Disentangling the relative contributions of shared ancestral polymorphism(i.e. incomplete lineage sorting) and recent gene flow to genetic differentiation is a major
challenge in population genetics research (Nielsen and Wakeley, 2001; Duthiel and Hobolth, 2012). In
short, low genetic differentiation may represent recent divergence, recent gene flow or both.
For Japanese stickleback, it is therefore possible that some freshwater populations wereoriginally established by the JA lineage but these have subsequently become isolated fromJA but not PA populations, leading to high and low differentiation respectively. We consider
this scenario an unlikely explanation for the genetic differentiation values estimated inthis study. For example, the higher proportion of private alleles across loci shared between

Fig. 4. Phylogenetic tree based on Nei’s DA (A) and DSW (B). Bootstrap values larger than 0.5 are shown at the branch. Solid circles, open circles,
and solid squares denote PA, JA, and freshwater populations, respectively.

Origin of Japanese freshwater stickleback populations

307

PA and freshwater populations suggests that the PA lineage is ancestral to freshwater
populations. Furthermore, differential amplification at the Stn76 locus also provides evidence of shared ancestry between freshwater and PA fish. Stn76 failed to amplify in all of
the JA populations analysed, but successfully amplified in all others. Based on the wholegenome sequences of an Akkeshi Japan Sea fish, the failure of amplification is likely due tothe existence of two nucleotide substitutions within the reverse primer sequence (K. Yoshida andJ. Kitano, unpublished data). As non-amplification of Stn76 is a shared characteristic common to
all JA fish in this study, such polymorphism probably occurred soon after the split between
JA and PA fish.
Nonetheless, a single freshwater population, Gifu, did group with the Japan Sea lineage
when K = 2 in the STRUCTURE analysis. However, there are several reasons we think the
Gifu population is likely derived from the PA lineage. First, genetic divergence between
the Gifu and JA populations is higher than between the Gifu and PA populations (Jost’s
D = 0.84 and 0.76 with Akkeshi JA and PA respectively). Furthermore, although in the tree
based on Nei’s DA, the Gifu population was located very close to the root of the Japanese
threespine stickleback, the tree based on DSW clearly demonstrated that the Gifu population is more genetically similar to the PA rather than to the JA populations. DSW generally
performs better than Nei’s DA for the resolution of the genetic relationships of distantly
related species (Goldstein et al., 1995; Terazaki and Nei, 1996). Because previous mitochondrial DNA
data indicate that the Gifu population likely diverged well before most contemporary
freshwater populations (Watanabe et al., 2003), a phylogenetic tree based on DSW might beappropriate to resolve the relationships between the Gifu population, the JA lineage, and
the PA lineage. In addition, the analysis of private alleles also supports the hypothesis that
the Gifu population is genetically similar to the PA populations. Furthermore, our recent
whole genome sequence data support the hypothesis that the Gifu population is geneticallymore similar to the PA than the JA populations (K. Yoshida and J. Kitano, unpublished data). Therefore, the Gifu population is also likely derived from the PA lineage.
Since none of the Japanese freshwater populations surveyed to date are derived from theJA lineage, this leads us to ask what factors might prevent freshwater colonization in thislineage. First, ninespine stickleback (genus Pungitius) inhabit freshwater lakes and rivers
connected to the Sea of Japan (Ikeda, 1933; Takahashi et al., 2001; Tsuruta and Goto, 2006). This contrasts
with the Pacific Ocean coast, where ninespine stickleback do not occur in freshwater
environments except on Hokkaido. Therefore, the ninespine stickleback may have occupied
these freshwater environments earlier than the threespine stickleback, precluding theircolonization. Second, JA fish may lack the key traits important for survival in freshwater environments. For example, previous studies demonstrated that JA fish kept infresh water had low survival rates in contrast to higher rates in PA fish (Honma, 1975; Yamada,2003). Differences in the survival rate in freshwater environments might be due to differencesin the ability of osmoregulation and fatty acid synthesis (A. Ishikawa and J. Kitano, unpublished data).
Understanding the factors that facilitate and constrain adaptive radiations is one of thecentral challenges in evolutionary ecology. Although novel niches are thought to promote
adaptive radiation, not all lineages are able to seize the ecological opportunity. Following
the last deglaciation, rising sea levels and shifting ice mass revealed multiple freshwater
environments. Several lineages of stickleback, smelts, and salmonids colonized these
newly formed lakes and rivers, resulting in impressive adaptive radiations (Bernatchez et al., 1999;
Taylor, 1999; Schluter, 2000; Kinnison and Hendry, 2004). In contrast, some lineages have failed to diversify.
The Japanese anadromous species pair provides a good demonstration of this. While308

Cassidy et al.

the PA lineage has repeatedly colonized freshwater environments, resulting in increased
ecological, morphological, and genetic diversity, the JA lineage has not radiated to any
substantial extent. Further comparative studies focusing on the PA and JA lineages will
provide important insight into the factors that facilitate and constrain adaptive radiation.
For example, further studies on the genetic and genomic mechanisms underlying thedivergence in freshwater adaptation between the PA and JA lineages and genomic
comparison of standing genetic variation between the two lineages will help to answer whysome lineages could colonize and utilize vacant niches to achieve adaptive radiation, whileothers could not.
ACKNOWLEDGEMENTS
We thank M. Kume, Y. Machida, and S. Takeyama for help with sampling, T. Kokita for providing
samples and helpful comments on the manuscript, A. Ishikawa and K. Yoshida for discussion, and the‘Share the Global Map’ Project for providing a map file. A.P. Hendry, C.L. Peichel, and an anonymous
reviewer made constructive comments on the manuscript. This research was supported by JST
PRESTO program, Grant-in-Aid for Scientific Research on Innovative Areas from the Ministry of
Education, Science, Sports, and Culture (23113007 and 23113001), and the Environment Research
and Technology Development Fund from the Ministry of the Environment (ZD–1203) to J.K. and
NIG Cooperative Research Program (2012–A63) to S.M. M.R. was supported by the Japanese Society
for the Promotion of Science on the Short-term Postdoctoral Fellowship for Foreign Researchers.
L.C. was supported by the NIG Intern Program 2012.

REFERENCES
Adachi, T., Ishikawa, A., Mori, S., Makino, W., Kume, M., Kawata, M. et al. 2012. Shifts inmorphology and diet of non-native sticklebacks introduced into Japanese crater lakes. Ecol.
Evol., 2: 1083–1098.
Avise, J.C. 2000. Phylogeography: The History and Formation of Species. Cambridge, MA: Harvard
University Press.
Ballard, L., Adams, P., Bao, Y., Bartley, D., Bintzler, D., Kasch, L. et al. 2002. Strategies forgenotyping: effectiveness of tailing primers to increase accuracy in short tandem repeat
determinations. J. Biomol. Tech., 13: 20–29.
Bell, M.A. and Foster, S.A. 1994. The Evolutionary Biology of the Threespine Stickleback. Oxford:
Oxford University Press.
Bennett, K.D. 2008. What do we mean by ‘refugia’? Quaternary Sci. Rev., 27: 2449–2455.
Bernatchez, L., Chouinard, A. and Lu, G. 1999. Integrating molecular genetics and ecology instudies of adaptive radiation: whitefish, Coregonus sp., as a case study. Biol. J. Linn. Soc., 68:
173–194.
Duthiel, J.Y. and Hobolth, A. 2012. Ancestral population genomics. In Evolutionary Genomics: Statistical and Computational Methods, Vol. 2: Methods in Molecular Biology (M. Anisimova, ed.),pp. 3–49. New York: Springer.
Earl, D.A. and vonHoldt, B.M. 2012. STRUCTURE HARVESTER: a website and program forvisualizing STRUCTURE output and implementing the Evanno method. Conserv. Genet.
Resources, 4: 359–361.
Evanno, G., Regnaut, S. and Goudet, J. 2005. Detecting the number of clusters of individuals using
the software STRUCTURE: a simulation study. Mol. Ecol., 14: 2611–2620.
Gillespie, R.G. and Baldwin, B.G. 2010. Island biogeography of remote archipelagoes. In The Theory
of Island Biogeography Revisited (J.B. Losos and R.E. Ricklefs, eds.), pp. 358–387. Princeton,
NJ: Princeton University Press.

Origin of Japanese freshwater stickleback populations

309

Glor, R.E. 2010. Phylogenetic insights on adaptive radiation. Annu. Rev. Ecol. Evol. Syst., 41:
251–270.
Goldstein, D.B., Ruiz Linares, A., Cavalli-Sforza, L.L. and Feldman, M.W. 1995. An evaluation ofgenetic distances for use with microsatellite loci. Genetics, 139: 463–471.
Hendry, A.P., Bolnick, D.I., Berner, D. and Peichel, C.L. 2009. Along the speciation continuum insticklebacks. J. Fish Biol., 75: 2000–2036.
Higuchi, M. and Goto, A. 1996. Genetic evidence supporting the existence of two distinct species in
the genus Gasterosteus around Japan. Envion. Biol. Fish., 47: 1–16.
Honma, Y. 1975. Excessively enlarged thyroid follicles of the threespine stickleback, Gasterosteus
aculeatus reared in freshwater. Jpn. J. Ichthyol., 21: 183–190.
Ikeda, K. 1933. Distribution and variation of sticklebacks. Zoolog. Mag. (Tokyo), 45: 141–173.
Jakobsson, M. and Rosenberg, N.A. 2007. CLUMPP: a cluster matching and permutation
program for dealing with label switching and multimodality in analysis of population structure.
Bioinformatics, 23: 1801–1806.
Jost, L. 2008. GST and its relatives do not measure differentiation. Mol. Ecol., 17: 4015–4026.
Kalinowski, S.T. 2004. Counting alleles with rarefaction: private alleles and hierarchical sampling
designs. Conserv. Genet., 5: 539–543.
Keenan, K. 2012. diveRsity: genetic diversity partition statistics and informative locus selection.
http:CRAN.R-project.org/package=diveRsity.
Kinnison, M.T. and Hendry, A.P. 2004. From macro- to micro-evolution: tempo and mode insalmonid evolution. In Evolution Illuminated (A.P. Hendry and S.C. Stearns, eds.), pp. 208–231.
Oxford: Oxford University Press.
Kitano, J., Mori, S. and Peichel, C.L. 2007. Phenotypic divergence and reproductive isolation
between sympatric forms of Japanese threespine sticklebacks. Biol. J. Linn. Soc., 91: 671–685.
Kitano, J., Bolnick, D.I., Beauchamp, D.A., Mazur, M.M., Mori, S., Nakano, T. et al. 2008a.
Reverse evolution of armor plates in the threespine stickleback. Curr. Biol., 18: 769–774.
Kitano, J., Mori, S. and Peichel, C.L. 2008b. Divergence of male courtship displays between
sympatric forms of anadromous threespine stickleback. Behaviour, 145: 443–461.
Kitano, J., Ross, J.A., Mori, S., Kume, M., Jones, F.C., Chan, Y.F. et al. 2009. A role for a neo-sex
chromosome in stickleback speciation. Nature, 461: 1079–1083.
Kume, M., Kitamura, T., Takahashi, H. and Goto, A. 2005. Distinct spawning migration patterns
in sympatric Japan Sea and Pacific Ocean forms of threespine stickleback Gasterosteus aculeatus.
Ichthyol. Res., 52: 189–193.
Kume, M., Kitano, J., Mori, S. and Shibuya, T. 2010. Ecological divergence and habitat isolation
between two migratory forms of Japanese threespine stickleback (Gasterosteus aculeatus). J. Evol.
Biol., 23: 1436–1446.
Losos, J.B. 2010. Adaptive radiation, ecological opportunity, and evolutionary determinism. Am.
Nat., 175: 623–639.
McKinnon, J.S. and Rundle, H.D. 2002. Speciation in nature: the threespine stickleback model
systems. Trends Ecol. Evol., 117: 480–488.
Meirmans, P.G. and Van Tienderen, P.H. 2004. GENOTYPE and GENODIVE: two programs forthe analysis of genetic diversity of asexual organisms. Mol. Ecol. Notes, 4: 792–794.
Mori, S. 1987. Geographical variations in freshwater populations of the three-spined stickleback,Gasterosteus aculeatus, in Japan. Jpn. J. Ichthyol., 34: 33–46.
Mori, S. 2003. Distribution of threespine stickleback in Japan: conserving stickleback populations.
In The Natural History of Sticklebacks (A. Goto and S. Mori, eds.), pp. 201–212. Sapporo:
Hokkaido University Press.
Nei, M., Tajima, F. and Tateno, Y. 1983. Accuracy of estimated phylogenetic trees from molecular
data. II. Gene frequency data. J. Mol. Evol., 19: 153–170.
Nielsen, R. and Wakeley, J. 2001. Distinguishing migration from isolation: a Markov Chain Monte
Carlo approach. Genetics, 158: 885–896.

310

Cassidy et al.

Ortí, G., Bell, M.A., Reimchen, T.E. and Meyer, A. 1994. Global survey of mitochondrial
DNA sequences in the threespine stickleback: evidence for recent migrations. Evolution, 48:
608–622.
Peichel, C.L., Nereng, K.S., Ohgi, K.A., Cole, B.L., Colosimo, P.F., Buerkle, C.A. et al. 2001. The
genetic architecture of divergence between threespine stickleback species. Nature, 414:
901–905.
Peichel, C.L., Ross, J.A., Matson, C.K., Dickson, M., Grimwood, J., Schmutz, J. et al. 2004. The
master sex-determination locus in threespine sticklebacks is on a nascent Y chromosome. Curr.
Biol., 14: 1416–1424.
Petit, R.J., Mousadik, A.E. and Pons, O. 1998. Identifying populations for conservation on the basisof genetic markers. Conserv. Genet., 12: 844–855.
Pritchard, J.K., Stephens, M. and Donnelly, P. 2000. Inference of population structure using multilocus genotype data. Genetics, 155: 945–959.
R Development Core Team. 2012. R: A Language and Environment for Statistical Computing.
Vienna, Austria: R Foundation for Statistical Computing.
Raymond, M. and Rousset, F. 1995. GENEPOP (version 1.2): population genetics software for exact
tests and ecumenicism. J. Hered., 86: 248–249.
Rosenberg, N.A. 2004. Distruct: a program for the graphical display of population structure. Mol.
Ecol. Notes, 4: 137–138.
Saitou, N. and Nei, M. 1987. The neighbor-joining method: a new method for reconstructing
phylogenetic trees. Mol. Biol. Evol., 4: 406–425.
Schluter, D. 1996. Ecological speciation in postglacial fishes. Phil. Trans. R. Soc. Lond. B, 351:
807–814.
Schluter, D. 2000. The Ecology of Adaptive Radiation. New York: Oxford University Press.
Shriver, M.D., Jin, L., Boerwinkle, E., Deka, R., Ferrell, R.E. and Chakraborty, R. 1995. A novel
measure of genetic distance for highly polymorphic tandem repeat loci. Mol. Biol. Evol., 13:
914–920.
Szpiech, Z.A., Jakobsson, M. and Rosenberg, N.A. 2008. ADZE: a rarefaction approach forcounting alleles private to combinations of populations. Bioinformatics, 24: 2498–2504.
Takahashi, H., Takata, K. and Goto, A. 2001. Phylogeography of lateral plate dimorphism in the
freshwater type of ninespine sticklebacks, genus Pungitius. Ichthyol. Res., 48: 143–154.
Takezaki, N., Nei, M. and Tamura, K. 2010. POPTREE2: software for constructing population trees
from allele frequency data and computing other population statistics with Windows interface.
Mol. Biol. Evol., 27: 747–752.
Taylor, E.B. 1999. Species pairs of north temperate freshwater fishes: evolution, taxonomy, andconservation. Rev. Fish Biol. Fish., 9: 299–334.
Terazaki, N. and Nei, M. 1996. Genetic distances and reconstruction of phylogenetic trees frommicrosatellite DNA. Genetics, 144: 389–399.
Tsuruta, T. and Goto, A. 2006. Fine scale genetic population sturucture of the freshwater andOmono types of nine-spined stickleback Pungitius pungitius (L.) within the Omono river system,
Japan. J. Fish Biol., 69 (suppl. B): 155–176.
Watanabe, K., Mori, S. and Nishida, M. 2003. Genetic relationships and origin of two geographic
groups of the freshwater threespine stickleback, ‘Hariyo’. Zoolog. Sci., 20: 265–274.
Waters, J.M. 2011. Competitive exclusion: phylogeography’s ‘elephant in the room’? Mol. Ecol.,20: 4388–4394.
Weir, B.S. and Cockerham, C.C. 1984. Estimating F-statistics for the analysis of population structure. Evolution, 38: 1358–1370.
Wootton, R.J. 1976. The Biology of Sticklebacks. London: Academic Press.
Wootton, R.J. 1984. A Functional Biology of Sticklebacks. London: Croom Helm.
Yamada, M. 2003. Evolutionary process and introgression of mitochondrial DNA in threespine
stickleback, Gasterosteus aculeatus, around Japan. PhD thesis, Hokkaido University.

Origin of Japanese freshwater stickleback populations

311

Yamada, M., Higuchi, M. and Goto, A. 2001. Extensive introgression of mitochondrial DNA foundbetween two genetically divergent forms of threespine stickleback, Gasterosteus aculeatus,
around Japan. Environ. Biol. Fish., 61: 269–284.
Yamada, M., Higuchi, M. and Goto, A. 2007. Long-term occurrence of hybrids between Japan Sea
and Pacific Ocean forms of threespine stickleback, Gasterosteus aculeatus, in Hokkaido Island,
Japan. Environ. Biol. Fish., 80: 435–443.

APPENDIX 1: Hardy-Weinberg equilibrium (HWE) test
Results of exact tests for HWE across loci in all sampled populations

Population

Stn170 Stn233 Stn64

Stn76 Stn159 Stn46

Stn90 Stn120 Stn278 Stn332 Stn384

AisakaAizu
Chimikeppu
Gensui
Gifu
Kussharo
Nishitappu
Shikotsu
TowadaAbashiri
Akkeshi JA
Akkeshi PA
ChoboshiHarutori JA
Harutori PA
HyotanSarufutsu
Watarichippu
Benkei
Notojima
Ogawarako
Onnechikappu

1.000
0.085
1.000
0.346
1.000
0.956
0.967
0.728
1.000
0.583
0.144
0.215
0.071
0.424
0.324
0.998
0.169
0.218
0.029
0.172
0.442
0.081

0.391
NP
0.017
0.680
0.200
0.863
0.520
1.000
1.000
NA
NA
0.550
0.883
NA
0.028
0.411
NA
1.000
NA
NA
NA
NA

0.306
NP
NP
1.000
NP
0.313
1.000
1.000
1.000
0.000*
0.007
0.773
0.713
0.000*
0.408
0.908
0.019
0.723
0.044
0.004*
0.000*
0.005

0.115
NP
0.747
1.000
1.000
0.668
0.325
0.217
0.221
0.738
0.676
0.814
0.559
0.169
0.190
0.748
0.755
0.189
0.335
0.099
0.190
0.044

0.862
0.021
0.001*
0.819
0.258
0.111
0.862
0.087
0.046
0.500
0.779
0.959
0.805
0.546
0.766
0.349
0.815
0.642
0.741
0.260
0.355
0.721

0.019
0.019
0.383
0.185
0.164
0.847
0.900
0.854
1.000
0.001*
0.008
0.301
0.215
0.625
0.274
0.287
0.408
0.855
1.000
0.016
0.373
0.455

0.479
NP
NP
1.000
NP
1.000
1.000
0.100
0.573
0.778
0.129
0.534
0.183
0.425
0.629
0.183
0.939
1.000
0.181
0.870
1.000
1.000

0.833
NP
1.000
1.000
NP
0.235
0.158
0.332
0.306
0.089
0.641
0.093
1.000
0.963
1.000
1.000
0.970
0.876
0.465
0.221
0.476
0.165

NP
NP
NP
NP
NP
0.490
1.000
NP
NP
0.278
0.548
1.000
1.000
0.436
1.000
1.000
0.826
0.753
1.000
0.587
0.467
0.023

0.307
1.000
1.000
1.000
0.711
0.816
0.241
1.000
NP
0.633
0.001*
1.000
0.254
0.014
0.034
0.152
0.414
0.544
0.788
0.014
1.000
0.148

0.170
NP
0.478
0.416
NP
0.219
1.000
0.641
1.000
0.479
1.000
0.984
0.740
1.000
0.273
0.338
1.000
0.625
0.122
1.000
0.729
0.183

Note: Bold text denotes loci significant after false discovery rate (B-Y method) has been accounted for; asterisk
denotes loci significant after Bonferroni correction (α = 0.05, αFDR = 0.0165, αBF = 0.005). NA, microsatellites were
not amplified (null allele). NP, not polymorphic.