The Migraine model most consistent with our data included two demographic changes, a change from an ancestral effective population size (Nanc) followed by an increase to a founder size (Nfound) until recently when the size stabilized at the present-day effective size (N). To interpret these sizes and times in the past, we will briefly review glacial history of the Barents Sea and nearby regions. Marine isotope stage (MIS) corresponds to a temporal index based on the ratio of the 16O and 18O isotopes in deep-sea core sediments; relatively higher levels of 18O reflect colder periods (e.g. Urey et al. 1951). From ~ 190–140 ka (MIS 6), the Late Saalian Glacial stage advanced far into northern Europe and western Russia (Hughes and Gibbard 2018), well past the limits of the Last Glacial Maximum (LGM) ~ 20–19 ka during the Late Weichselian (MIS 2; 25–15 ka) (Svendsen et al. 2004). There is evidence that much of the land-locked Arctic Ocean was covered with an ice sheet at MIS 6 (~ 140 ka) and high points on the ocean floor probably pinned the > 1000 m thick sheet (Jakobsson et al. 2016). This means that there would have been little macrofauna at that time in the western Arctic Ocean and colonization of species like navaga must have occurred subsequently.
Following the Late Saalian period was an interglacial, the Eemian period (MIS 5e; ~130–115 ka), which may have been warmer than our present climate (Sakari Salonen et al. 2018). Next came the Weichselian period, during which glaciers advanced three times: (1) the Early Weichselian (100–80 ka; MIS 5d), (2) the Middle Weichselian (60–50 ka; MIS 3) and (3) the Late Weichselian (25–15 ka; MIS 2). Each was followed by deglaciation (Svendsen et al. 2004). The Early Weichselian glaciation was the most extensive; the Kara, Barents, and Scandinavian ice sheets coalesced to extend from the Scandinavian Peninsula to the Taymyr Peninsula, which is the western boundary of the Laptev Sea. Grounded ice completely covered the Barents and Kara Seas from the coastline to the continental shelf. Inside the Russian ice sheet, large proglacial lakes formed, which probably drained through the Volga River to the Aral Sea (Soulet et al. 2013). When the glaciers receded, an enormous volume of water erupted through the southern Barents Sea and through the Baltic Sea (Mangerud et al. 2004). A combination of the grounded ice and huge water flow through the shallow coastal waters precludes the existence of navaga at our sampling sites during this period. During the Middle Weichselian, grounded ice sheets again were continuous from the Scandinavian Peninsula to the Taymyr Peninsula and proglacial lakes formed south of the Russian glaciers. At the end of this period, the proglacial lake waters probably flowed north. During the Late Weichselian, coalesced ice sheets extended through the Barents and Kara seas, but there appeared to be only peripheral proglacial lake formation (Mangerud et al. 2004). The regional ice volumes varied during the three Weichselian glaciations. Ice was thickest in the east during the Early Weichselian; whereas, there was more extensive glaciation from the Scandinavian Ice sheet during the Late Weichselian (Svendsen et al. 2004). Nevertheless, there is evidence of grounded ice throughout the Barents Sea (Landvik et al. 1998) during all three periods.
The demographic estimates from the microsatellite data can be interpreted in terms of the glacial history. During the Eemian interglacial, which followed the Late Saalian when much of the Arctic Ocean was ice covered, nearshore habitat became available. Our analysis detected a first demographic event that was ~ 140 ka, which is consistent with the Eemian interglacial. The effective population size at that time was about 2700 individuals. The ancestral population preceding that time was much smaller than that of the present day (Fig. 2, Table 2). The results of the model we applied in Migraine indicated that the current population size (about 26,000) had been attained by the time of the second demographic event. That means that the population had expanded about 10-fold since the Eemian. Constrained by the models available in Migraine, the population growth was modeled assuming exponential growth. However, it is unlikely continuous growth occurred because appropriate habitat did not become instantaneously available, rather it probably expanded slowly. Habitat in the Barents Sea was unavailable until ice retreated starting ~ 15 ka and the southern connection between the Barents and Kara seas was not available until ~ 10 ka (Fig. 3). Also, most of the rise in sea level was not complete until ~ 7 ka (Lambeck et al. 2014). It would be expected that habitat did not actually become available to navaga until the ecosystem expanded to include an appropriate complex of prey populations, which would lag behind simple access to the shallow nearshore marine environment. From the known present-day distribution, it seems likely that the colonization occurred as a westward movement from the Laptev or East Siberian seas. Because the western Russian Arctic coast was not continually ice-covered during the Weichselian advances, that area may have been a refugium. The estimate of the present-day effective population size from the coalescent model of Migraine was ~ 26,000, which was included in the range (1778–33,952) estimated from the linkage disequilibrium method. In contrast, the census size exceeds 4 million fish. The estimate from Migraine makes use of the distribution and suite of alleles to construct an historic estimate; whereas, NeEstimator looks at the linkage disequilibrium for pairs of alleles, but does not use the information provided by the sizes of the alleles. Regardless, the census size far exceeds effective numbers, possibly as a consequence of the recent colonization compounded by large variance in family size and inter-annual variation in the number of recruits.
Our last question is the source of the colonization of navaga following the LGM. Nucleotide sequence data from the mitochondrial cytochrome B gene of navaga and its sister taxon, saffron cod (E. gracilis), suggest that they diverged ~ 2.32 Ma (million years ago), as did several other vicariant species pairs (Laakkonen et al. 2021). The Pliocene cooling that took place 3.5–3.1 Ma (e.g. De Schepper et al 2015) was followed by a large glacial advance between ~ 3.6 and 2.4 Ma (Matthiessen et al 2009), after which cooling, relative to the earlier Pliocene and Miocene, continued through the Pleistocene to the present. Prior to complete closure of the Isthmus of Panama, the connection between the Pacific and Atlantic oceans had a depth of at least 1200 m (O’Dea et al. 2016). Subsequently, global oceanic current patterns were radically altered (e.g. De Schepper et al. 2015). One of the consequences of the changes in oceanic currents was the reversal of flow through the Arctic from the Pacific Ocean to the Atlantic Ocean during interglacial periods, when the water level was high enough to inundate the Bering Strait (e.g. De Schepper et al. 2015).
The mitochondrial nucleotide sequence data suggest that there was no contact between navaga and saffron cod subsequent to their divergence (Laakkonen et al. 2015). Biogeographical patterns deduced from mitochondrial gene sequences suggest that 34 taxa, which include mollusks, crustaceans, and fish also have a deep vicariance between the Pacific and Atlantic oceans (Laakkonen et al. 2015). Although E. nawaga and E. gracilis have slightly diverged morphologically; they can be difficult to distinguish without close scrutiny that may include dissection (Vasil’eva 1997).
During the Pleistocene, multiple glacial advances blocked access through the Bering Strait (Miller et al. 2020). As in the Late Pleistocene (the Saalian and Weichselian glacial cycles), the Arctic Ocean of northern Europe and northwestern Russia was completely covered by glaciers, as was the coastline of Greenland and North America. The Arctic coastal areas of eastern Siberia, which included much of the large continental shelf, were not glaciated (e.g. Batchelor et al. 2019), but were probably covered by seasonal ice. The resulting large coastal plain was grasslands, which supported grazers such as wooly mammoths and horses (Sher et al. 2005). Also during the late Pliocene and early Pleistocene, because of uplift in eastern Asia, Siberian rivers (Lena, Ob, and Yenisey) began to drain into the Arctic Ocean (e.g. Ma et al. 2021). Those rivers now contribute much of the total freshwater flow into the Arctic Ocean. The estuaries that were created could have served as favorable navaga habitat. It seems likely that the Laptev and Eastern Siberian seas may have provided refuge for navaga during glacial advances that occurred subsequent to the separation of navaga and saffron cod lineages.
What do we know about the ranges of navaga and saffron cod at present? There appears to be a several thousand kilometer disjunction, from the western Laptev Sea to the Chukchi Sea. The documented eastern limit for navaga is the Khatanga Bay in the western Laptev Sea (Ulchenko et al. 2016). Records of more eastern navaga are erroneous or undocumented. The most cited original report (Borisov 1928) of navaga for the Laptev Sea (the Lena River delta) was a misidentification of an Arctogadus borisovi specimen (Svetovidov 1948). There are no recent records of navaga in field studies in the Lena drainage, the Lena Polynya, the Laptev Sea shelf, or in the East-Siberian Sea (Gukov 1999; Chernova 2015; Orlov et al. 2020a, b).
The western limit of saffron cod is Kolyuchinskaya Bay in the western Chukchi Sea (Rendahl 1931). Saffron cod were not listed for Chaunskaya Bay in the western East Siberian Sea (Neelov 2008; Chernova 2022) nor caught in trawl surveys on shelf of the East Siberian Sea (e.g. Glebov et al. 2016; Orlov et al. 2020c). Although saffron cod are included in a few compilations of fish species of the East Siberian and Laptev seas), there is no documentation. Thus, at present navaga and saffron cod seem to be absent in the area between the Khatanga Bay (western Laptev Sea) and Kolyuchinskaya Bay (western Chukchi Sea).
In the coastal areas of the Siberian Arctic, the modern climate is colder than it was in the late Pleistocene – early Holocene. This process is well documented by the succession of paleoflora (Svitoch 1980). During the climate optimum, the coastal region was occupied by forest and shrub tundra with thickets of birch, alder, and dwarf forests. As a result of a climate cooling that started 5–4 ka, this vegetation was replaced by arctic-type tundra. Paleoichthyology provides similar evidence. Fossil remains of freshwater fish were identified in sediments (~ 12 ka) in the region of Chaunskaya Bay (Nazarkin 1992). Several species (e.g. Crucian carp Carassius sp., lake minnow Phoxinus percnurus (= Rhynchocypris percnura), longnose sucker Catostomus catostomus, and European perch Perca fluviatilis) no longer inhabit the area, presumably because of Holocene cooling (Nazarkin 1992).
Perhaps Holocene climate cooling explains the modern absence of the navaga and saffron cod in the coastal Siberian water. The Laptev and East Siberian seas are among the harshest Arctic areas. The duration of the ice-free period ranges from 1–1.5 months in the northeast to 2.5–3 months in the southeast. The greatest influx of solar radiation falls on May-June, when the sea is still covered with ice and 60% of the incoming radiation reflects into the atmosphere (Ecological Atlas 2017). Warm currents of neither the Pacific nor Atlantic oceans penetrate into Siberian shelf areas (Anderson and Macdonald 2015). Although water temperatures have increased during the last 20 years (Golubeva et al. 2021), the East Siberia shelf is still less saline (< 20‰) and colder (1–5 C) in the summer than the Kara, Barents, and White sea coasts and freshwater intrudes down to ~ 25 m along much of the coast. The enormous rivers also deliver substantial amounts of clay and silt, which characterize the very broad shelf of the East Siberian Sea (e.g. Kokarev et al. 2021). The result of the physical environment is an ecosystem with low productivity. All these factors make habitat unfavorable for both navaga and saffron cod. This cold Siberian “spot” may limit the distribution of saffron cod from the Pacific along the Asian coast. It also may have caused the disappearance of navaga from this area. Navaga may have been distributed more widely during the boreal optimum, eastward of the Kara-Barents seas, but after the Holocene cooling, it could have been concentrated along the western Siberian coast; the modern range of navaga may be a recent refuge.