Comparison among markers and demographic interpretations
The extreme difficulties of sampling in Antarctica resulting in relatively low sample sizes and coverage, prevent any detailed conclusions about Trematominae population structure. Some aspects deserve discussion, however, especially due to the potential effect on demographic inferences. Whereas analysis of the S7 locus detected no deviations from panmixia in both species, analysis of mtDNA suggested a distinctive South Shetland population of T. bernacchi, which was previously postulated . Such discrepancies between nuclear and cytoplasmic markers have been more or less frequently encountered in other studies, and may result from a number of mechanisms. Whereas some studies clearly stand in support of deterministic explanations, such as selective pressures or different evolution of sexes [66, 67], other cases (including Antarctic toothfish ) were explained by a simpler stochastic effect of different coalescence dynamics of mtDNA and nuclear markers [e.g. [69–73]]. A four-fold smaller effective population size of mtDNA relative to nuclear loci makes it more sensitive to population oscillations [ and citations therein], and suggests that genetic drift will result in a faster population subdivision relative to a nuclear locus.
We prefer the latter explanation given the lack of evidence of sex-specific dispersal in Trematominae and no apparent signs of either direct or background selection on polymorphism in either locus. Trematomus bernacchi, in common with other benthic fish  shows a pronounced population structure on a continent-wide scale. Grounded ice sheets did not cover the whole shelf in the Ross Sea during LGM [75, 76], which is therefore considered as a potential refugium for benthic biota in particular . Together with a decrease of oceanic current activities during glacial periods  reducing passive larval dispersal [79, 80], our data are concordant with fragmentation into separate refugia with subsequent independent population expansions. Our sampling is far from complete and the following interpretations are only very approximate. The observed patterns and assumption of geographic homogeneity inherent with the demographic inference tools, at least justify the analysis of the Ross Sea population as a separate entity. The South Shetlands population was not analysed separately due to the low sample size from this location.
Despite generally concordant signals of population expansion from both loci, demographic analyses indicated another potential conflict between mtDNA and nuclear markers: 1) We couldn't significantly reject the constant-size model in traditional methods of Tajima's D and Fu's Fs in the S7 datasets despite negative values, 2) estimates of g were lower in nuclear genes and 3) inferred expansion dates were much more ancient in S7. Locus-specific selective sweeps may explain such discrepancies [ and citations therein], but our results make such an explanation unlikely. Whereas Tajima's D and Fu's Fs are sensitive to demographic oscillations, we detected no signs of selection in the McDonald-Kreitman or HKA tests, which are not expected to be affected by population size changes. Furthermore, these traditional measures seem more conservative than ML methods, which also use phylogenetic information [19, 61]. Our interpretation of concordant expansion signals is further strengthened by the similar results of Galtier et al. test for both datasets of benthic species. In addition, significantly positive population growth estimated independently from both loci was supported by the results of a combined analysis, which is known to effectively reduce the upward bias in estimates of g .
In our opinion, contrasting signals from cyt b and S7 markers with respect to expansion rates, times and the significance of MD, Tajima's D and Fu's Fs, are better explained by the different properties of each marker and population histories that were not identical to the simple models assumed by demographic methods. Given that g is scaled by the locus-specific mutation rate , lower rates in S7 likely come from an approximately four-fold difference in mutation rate of both markers (M. Kuhner and L. Smith, pers. comm.). It also suggests that a faster-evolving locus like cyt b is likely to incorporate traces of more recent demographic events . Furthermore, demographic histories of the studied species were probably more complex than simple exponential or stepwise expansion due to repeated climatic cycles during the Pleistocene. Rogers and Harpending  and Rogers  showed that in cases of complex population histories, the parameter τ correlates with the initial expansion, which would obscure the effects of later events (including bottlenecks) for some time. Nuclear markers with longer coalescent times relative to mtDNA are expected to accumulate the information of demographic history over longer periods and indicate older events. It is also possible that, compared to the cyt b, temporal changes in population growth rate attenuated "more ancient" signals of population expansion in S7 to a greater extent.
Clear discordance between cyt b and S7 was nonetheless observed in the T. newnesi data, where positive values of Tajima's D and Fu's Fs (significant Tajima's D) may suggest population subdivision (albeit not detected in AMOVA, Φst and HW analyses) or balancing selection in S7 (but see the discussion above). We also reanalysed T. newnesi data with exclusion of the South Shetland population to avoid the effect of any undetected population structure, but obtained almost identical results (data not shown). Alternatively, positive values may also indicate a drop in population size, which is further supported by the method of Wakeley and Hey and insignificance of g estimates. Although combined analysis of both loci suggested positive growth in this species, it was the lowest among all species and our data probably do not allow us to make strong conclusions about the demographic history in this case.
When looking at the time-scale of faster-evolving mtDNA, which also allows the comparison with other published studies, the expansion events in pelagic species seem more ancient compared to benthic ones (the expansion of T. newnesi – if any – has the oldest date) although confidence intervals were quite large and partly overlapping especially between the P. borchgrevinki and T. bernacchi complete dataset. Unresolved estimates of τ for T. pennelli and the T. bernacchi Ross Sea population based on Schneider and Excoffier's  method (compare with lower τ values estimated from whole T. bernacchi data where the MD was skewed more to the right) may suggest that there were too few informative sites in our data for this method to provide reliable estimates. In any case, we have shown that lower C.I. for expansion event in P. borchgrevinki is older than the age of whole genealogy of the Ross Sea population of T. bernacchi (the age of T. pennelli is even younger). Since all inferences of demographic histories were based in one or the other way on mutation spectre in observed data, it is not possible, in our opinion, to detect expansion times older than the data themselves, supporting a more recent expansion of benthic species.
In summary, while the demographic history of T. newnesi remains somewhat obscure, the data from remaining three species strongly indicate that they underwent significant population growth events during the Pleistocene. Given the results of Galtier et al. test, the proposed population expansions may be linked to recovery from recent strong bottlenecks at least in case both benthic feeders (this method does not indicate the beginning of the event however, as the same bottleneck effect may come from short and strong population reduction or a long but weak one). Given the time estimates from both loci, we assume that populations of these species were fluctuating during much of the Pleistocene. The more recent expansion times of benthic feeders assessed by means of mtDNA may suggest different reactions of each group to historical climatic events.
Using the comparative approach , there is no apparent distinction between benthic and pelagic species with respect to the expansion rate, since ML values of g were overlapping between both groups. Thus, we may not confirm the hypothesis that population oscillations were more dramatic in benthic species. Our data suggest that the benthic T. bernacchi stands at one extreme of the expansion rate estimates, whereas the pelagic T. newnesi stands at the other.
Comparisons to other Southern Ocean organisms
Increasing numbers of studies of Northern Hemisphere marine fish suggest their populations have generally undergone important population size oscillations, which are interpreted in relation to Pleistocene changes of currents, sea-level or temperature [e.g. [83–89]]. Some species-specific ecological characteristics may nonetheless buffer the impact of climatic fluctuations on populations of marine organisms, while others, such as restriction to habitats prone to dessication or greater specialisation in diet and habitat preferences, make them more vulnerable to such changes [e.g. [90–93]].
In constrast, the scarcity of population genetic and phylogeographic studies on Southern Ocean marine organisms limits our understanding of Southern populations. Our finding of more recent expansion dates in benthic Trematomus suggest that habitat preference may have modified the impact of Pleistocene climatic shifts in agreement with the glacial-mediated habitat eradication hypothesis. However, absolute dating of demographic events and their relation to particular climatic events may be compromised by errors associated with molecular clock calibration and should be considered as approximate. On the other hand, such estimates may provide useful comparisons among species. Provided our calibration is correct, our data seem to provide a surprisingly good match of T. pennelli expansion to the onset of the last deglaciation, beginning about 25 kya but possibly as early as 30 kya in some parts of Antarctica . The T. bernacchi expansion also probably happened during the last major glacial period (or during the last ~50 kya when considering the 3rd cyt b position dataset). The best estimates from our data suggest that major population expansions of both pelagic Trematominae species by far predated the Holocene, and even the last interglacial, although, strictly speaking, the expansion of P. borchgrevinki might have occurred as recently as 70 kya according to 95% C.I..
A similar picture is suggested by data published so far, since other pelagic organisms, namely the krill Euphausia superba  and the silver fish Pleuragramma antarcticum , have most likely experienced populations expansions and/or bottlenecks by far predating the Holocene, dated to between 205 – 304 kya according to alternative molecular clocks and around 120 kya (48 – 142 kya according to 95% C.I.) in each species respectively. In their study of three benthic icefish (Chionodraco), Patarnello et al.  found high haplotype diversities, low divergences and negative values of Tajimas' D otherwise supportive for population bottlenecks. However, the insignificant values found in this study prevented the authors from reaching conclusions about the demographic history of these icefish. Reanalysing these data, we found in all three species smooth unimodal mismatch distributions and significantly negative values of Fu's Fs in C. hamatus. In the remaining species, both Tajima's D and Fu's Fs were significant at 10% level. We also found significantly positive values of g on the order of 102 in all three species (data available upon request). Using the same substitution rate as in , the population expansions in two of the three benthic species seem more recent compared to previously published data on pelagic feeders – most likely values suggested their onset approximately at 90, 47 and 25 kya in C. hamatus, C. myersi and C. rastrospinosus, respectively (95% C.I. suggested upper dates at 180, 95 and 90 kya, respectively).
To compare these data with our results we have to synchronise molecular clocks among the above studies [16, 74] and this paper. This is because we calibrated the clocks by the separation of Channichthyidae and Nototheniidae dated to 24 Mya according to the fossil-calibrated molecular clock , while Zane et al.  used as calibration the split between Eleginops maclovinus and the rest of Notothenioidei dated approximately to the same time (22–25 Mya) because supposedly corresponding to the formation of the APF (surprisingly, that divergence was set to about 40 Mya using the fossil calibration ). On the one hand, the calibration in Zane et al.  refers to previous reports [5, 95] using a calibration point of "perciform radiation" at 60 Mya, which is highly tentative, since we actually don't know what a perciform is . This consequently led to a correspondence between the formation of the APF and the separation of Eleginops from the rest of the notothenioids. In our opinion, this event rather correlates with the diversification of AFGP-bearing notothenioids (i.e. the divergence time between Pleuragramma and Chionodraco), which corresponds to the data on Notothenioidei [26, 97]. On the other hand, Near  based the calibration of his molecular clock on a single point with a fossil record and there have been doubts that this fossil really was a notothenioid. Both calibrations therefore suffer from some uncertainties.
Applying either one or the other calibration point to all fish datasets, we generally observed more recent expansion events in the benthic species relative to the pelagic organisms, except for C. hamatus. This indicates that habitat preferences affect the response to climatic shifts. Using the calibration date as in Zane et al.  (resulting in 1.6 to 1.8 times faster mutation rates compared to our calibration), the data would further suggest that population expansions of both benthic Trematomus correlate with the last glacial retreat.
On the other hand, the non-concordant signals in the rate and onset of population expansions within both groups of species complicate a clear-cut general conclusion. In addition to problems in using molecular clock estimates to compare the results for organisms as unrelated as fish and krill, it is possible that the sensitivity of species to glacial oscillations is determined by more complex traits than simple distinction between benthic and pelagic lifestyles. For example, some benthic species, like T. bernacchi, prefer shallow waters, feed mostly on benthic prey and use large sponges as spawning grounds, and thus require undisturbed habitats , while other benthic species such as Chionodraco spawn on rocks, consume large amounts of Euphausiids and pelagic fish, prefer greater depths and therefore might not necessarily be so vulnerable to glacial expansions.