INTRODUCTION
The use of parasite species as biological tags (BTs) is one of several approaches applied to identify fish stocks and to trace movements of fish populations (MacKenzie and Abaunza, Reference MacKenzie, Abaunza, Cadrin, Friedland and Waldman2005). Fish stock identification is essential for fishery sciences since it is a prerequisite for the assessment of exploited resources and their population dynamics, in order to provide information necessary for their proper management (Cadrin et al. Reference Cadrin, Friedland, Waldman, Cadrin, Friedland and Waldman2005). Over recent decades, the use of parasites as BTs has been increasingly applied in such sense, being at present a well-established methodology (Timi and MacKenzie, Reference Timi and MacKenzie2015).
Parasite species that are ideal as BTs are characterized by showing significantly different levels of infection among localities, persisting for long periods in the host, being easily detected and identified, and not causing serious damage to the host (MacKenzie and Abaunza, Reference MacKenzie, Abaunza, Cadrin, Friedland and Waldman2013; Catalano et al. Reference Catalano, Whittington, Donnellan and Gillanders2014). Some myxozoans meet the requirements, being host/site specific and easily identified by their distinctive spore shapes, with some species known to persist for years in their hosts (Lom and Dyková, Reference Lom and Dyková1992). Indeed, many myxozoan species have been successfully used as tags for stock identification of valuable fish species around the world (Moles and Short Reference Moles and Short1989; Khan and Tuck, Reference Khan and Tuck1995; Urawa and Nagasawa, Reference Urawa and Nagasawa1995; Larsen et al. Reference Larsen, Hemmingsen, MacKenzie and Lysne1997; Campbell, Reference Campbell2005; MacKenzie et al. Reference MacKenzie, Kalavati, Gaard and Hemmingsen2005). Moreover, this group of parasites proved to be a very effective BT on the coast and the open seas of the Pacific Ocean in North America (Marcogliese and Jacobson, Reference Marcogliese and Jacobson2015). The use of parasites as BTs has been successfully applied to delineate fish stocks in the South West Atlantic (Timi, Reference Timi2007; Cantatore and Timi, Reference Cantatore and Timi2015), only few studies using myxozoan species as BTs have been carried out in this region, with most of them corresponding to stock assessments of deep water fishes inhabiting the southern boundaries of the continental shelf (Agnew et al. Reference Agnew, Marlow, Lorenzen, Pompert, Wakeford and Tingley2003; George-Nascimento et al. Reference George-Nascimento, Moscoso, Niklitschek and González2011; MacKenzie et al. Reference MacKenzie, Brickle, Hemmingsen and George-Nascimento2013). Nonetheless, due to the oceanographic complexity of Southwest Atlantic at lower latitudes, it is expected that the discriminatory capability of myxozoans is at least equivalent to those demonstrated for several helminth species in that region (Cantatore and Timi, Reference Cantatore and Timi2015).
In this region the Argentine hake, Merluccius hubbsi Marini, 1933 is one of the principal sources of income in the Argentine and Uruguay fisheries (Bezzi et al. Reference Bezzi, Renzi, Irusta, Santos, Tringali, Ehrlich, Veis, Sánchez and Bezzi2004). This temperate-cold species is widely distributed in the Southwest Atlantic Ocean, inhabiting Brazilian, Uruguayan and Argentine waters, from 22° to 55°S at depths ranging between 50 and 500 m (Cousseau and Perrotta, Reference Cousseau and Perrotta2004).
Based upon evidence from ichthyological studies, two different stocks of M. hubbsi, limited by the 41°S parallel, have been identified in Argentine waters and considered for management purposes: the northern stock called Bonaerense and the southern one, known as Patagonian (Bezzi et al. Reference Bezzi, Verazay, Dato, Alheit and Pitcher1995). Each stock has its main spawning and nursery area and season; the northern area situated in the Argentine–Uruguayan Common Fishing Zone at the Rio de la Plata estuarine front during the austral autumn–winter (May–July) and the southern area in coastal waters of Chubut Province during the austral spring–summer (with a main peak in December and January) (Macchi et al. Reference Macchi, Pájaro and Ehrlich2004, Reference Macchi, Pájaro and Dato2007). It has also been established that hakes of both stocks perform trophic and reproductive seasonal migrations between the coast (near the 50 m isobath) and deeper waters (more than 100 m depth) (Pájaro et al. Reference Pájaro, Macchi and Ehrlich2005; Macchi et al. Reference Macchi, Pájaro and Dato2007). The existence of an independent third group, restricted to the San Matías Gulf basin (41°–42°S, 63°–65°W), has also been proposed based on fishery biology and population studies (Ehrlich, Reference Ehrlich1998; Sardella and Timi, Reference Sardella and Timi2004; Machado-Schiaffino et al. Reference Machado-Schiaffino, Juanes and Garcia-Vazquez2011).
Two previous studies have used parasite tags to discriminate stocks of M. hubbsi in the region; the first by MacKenzie and Longshaw (Reference MacKenzie and Longshaw1995) studied the potential value of parasites as BTs for Argentine hake in southern deep waters of Patagonia. However, as stated by the authors, their samples consisted of small numbers of fish of mixed sizes and ages, which precluded them from reaching reliable conclusions regarding stock identification and migrations of this host species. Nevertheless a number of parasite species were pointed out as the most promising BTs; among them, the authors identified the myxozoans Alatospora merluccii Kalavati, Longshaw and MacKenzie, 1995, Kudoa rosenbuschi (Gelormini, 1944), Myxidium baueri Kovaleva and Gaevskaya, 1982 and Myxoproteus meridionalis Evdokimova, Reference Evdokimova1977, the four myxozoan species so far known to infect the hake. The second by Sardella and Timi (Reference Sardella and Timi2004) identified three stocks of hake in the Argentine Sea using parasites as BTs, one composed of hakes from the Argentine–Uruguayan Common Fishing Zone, a second stock in San Matías Gulf and a third comprising samples from both San Jorge Gulf and the Patagonian Shelf. However, of the two myxozoan species found in this study, only K. rosenbuschi was useful as a discriminating species. The other species found, M. meridionalis, was not included in analyses because not all hakes were examined for it (Sardella and Timi, Reference Sardella and Timi2004). Based on these findings, myxozoans could be a good ancillary tool for delineating hake stocks, even at smaller spatial scales than in previous studies.
The main objectives of the study were to identify myxozoan species that parasitize M. hubbsi in the northern and central Argentine Sea and to assess the feasibility of using them as BTs to (i) discriminate among different hake stocks, and (ii) to trace movements of hake populations in the region.
MATERIALS AND METHODS
Fish samples and parasites inventories
In total, 213 specimens of Argentine hake were caught in two regions: in Patagonian waters (Chubut Province, 42°–46°S/63°–66°W) during a research cruise in January 2009 and by commercial trawlers in January 2014, and off the coast of Buenos Aires Province (36°–39°′S/53°–57°W) by commercial trawlers in July 2009 and May 2014 (Table 1).The two most recent samples, caught 5 years later, were included in order to evaluate temporal stability of parasite prevalence. As M. hubbsi displays different migratory patterns at both sampling regions concentrating in two separate and well-known reproductive areas (Bezzi et al. Reference Bezzi, Verazay, Dato, Alheit and Pitcher1995), the samples were collected during the reproductive season to ensure their belonging to each individual stock. Moreover, data on the exact geographical coordinates of capture available only for hakes caught during the research cruise (January 2009), allows comparative analyses among samples Patagonian 1, 2 and 3 to test for differences in parasitefaunas obtained in the regions of origin previous to reproductive migrations.
Composition of samples of Merluccius hubbsi Marini, 1933 (Gadiformes: Merlucciidae) in two regions of the South West Atlantic

Fish were kept frozen at −18° C until examination. Each specimen was measured (total length), sexed and thoroughly examined for myxozoan infections. Any tissue lesion or cyst observed under a stereoscopic microscope, as well as smears of brain, liver, kidney and spleen, and contents of gall and urinary bladders were examined with the aid of a light microscope equipped with differential interference contrast. Myxozoans were collected and preserved in 10% formaldehyde solution for subsequent morphological identification following the recommendations made by Lom and Arthur (Reference Lom and Arthur1989) and Lom and Dyková (Reference Lom and Dyková1992).
Statistical analyses
Prevalence was calculated according to Bush et al. (Reference Bush, Lafferty, Lotz and Shostak1997) for each myxozoan species at each sample. Sterne's exact 95% confidence limits were calculated for prevalence using Quantitative Parasitology 3·0 software (QP3·0) (Reiczigel, Reference Reiczigel2003; Reiczigel and Rózsa, Reference Reiczigel and Rózsa2005).
The effects of fish total length (L, in cm), year (Y, 2009 and 2014), region of origin (R, Bonaerense -BA- and Patagonian -PAT-) and their interactions on the occurrence of each myxozoan species were evaluated using generalized linear models (GLMs) (McCullagh and Nelder, Reference McCullagh and Nelder1989). Samples Patagonian 1, 2 and 3 (see Table 1), obtained from a research vessel as three different fishing hauls from the same region and year, were pooled in order to give the same treatment to all samples in the analysis. In the case of Fabespora, which was absent from Bonaerense samples (see results, Table 2), region of origin was excluded as a predictor variable because the complete separation of samples as a result of a ‘too good predictor’ produces extreme standard errors (s.e.) that are indicative of poor convergence of a model, resulting in unreliable inferences (Bolker et al. Reference Bolker, Brooks, Clark, Geange, Poulsen, Stevens and White2009).
Prevalence % (P) of myxozoan species found in Merluccius hubbsi Marini 1933

CI, Confidence interval.
Data on the exact geographical coordinates of capture were available only for hakes caught during the research cruise (January 2009), so analyses at a small spatial scale were only possible for these samples. GLMs were repeated excluding the variable year and replacing region of origin by haul position (Patagonian 1, 2 and 3) to evaluate its effects on the presence of each myxozoan species.
The variable sex was not included in all former analyses because not all combinations of year (2009 or 2014) and region of origin (Bonaerense or Patagonian) had equivalent female:male ratios (Table 1), resulting therefore in collinearity between sex and year/region of origin combinations. Nonetheless, the characteristics of the sample Patagonian 4 in particular (sample size, comparable sex proportions and parasite prevalence values) allowed the evaluation of the effect of sex on the occurrence of each myxozoan species within this sample. Therefore, new GLMs, including length, sex and their interactions as predictor variables were fitted and analysed (see below).
Following Mazerolle (Reference Mazerolle2006), for each myxozoan species in every analysis, a set of alternative candidate models was selected among all possible models by choosing those combinations of predictor variables and their interactions considered a priori as biologically meaningful (see details of candidate models in Tables 3 and 5). Also, null models (those that include no predictor variables) were fitted to assess if any of the input variables has an effect on the parasite occurrence. The presence–absence of the parasite species was the response variable and all models in each set had a binomial error distribution and a logit link function. The fit of each subglobal model, defined as the most complex model in each set, was assessed visually by building logistic regression quantile–quantile plots and partial residual plots for continuous covariate (length) (Landwehr et al. Reference Landwehr, Pregibon and Shoemaker1984; Zuur et al. Reference Zuur, Ieno, Walker, Saveliev and Smith2009).
Alternative candidate GLMs ranked according to AICc to evaluate the effects of predictor variables (L total length, Y year and R region of origin) on the occurrence of three myxozoan parasites (Kr – Kudoa rosenbuschi, Mm – Myxoproteus meridionalis, Fsp – Fabespora sp.) of Merluccius hubbsi in two regions of the South West Atlantic

ID, model identification; K, number of explanatory parameters in the model; AICc, second-order Akaike Information Criterion; Δi, AICc(best)–AICc(min); wi, Akaike weight; acc wi, cumulative Akaike weight; LLlog-lokelihood; w′i recomputed wi for model averaging; over the dotted line models with acc wi < 0·9; in bold non nested models used for multimodel inference.
a Subglobal model of each set.
b Null model of each set.
The best model of each set was selected using the Information Theoretic (IT) model comparison method (Anderson et al. Reference Anderson, Burnham and Thompson2000; Burnham and Anderson, Reference Burnham and Anderson2002), which provides a framework that allows a simultaneous assessment of several, a priori competing biological hypotheses (Burnham and Anderson, Reference Burnham and Anderson2002; Johnson and Omland, Reference Johnson and Omland2004). The process of identifying the best model involves trading off model bias against model precision and IT achieves this by using appropriately constructed criteria to compare models (i.e. Akaike Information Criterion) (Burnham and Anderson, Reference Burnham and Anderson2002).The relative fit of each candidate hypothesis was quantified by mean of second-order Akaike Information Criterion (AICc) which is adjusted for small samples. The AICc was used because it approximates AIC at large sample sizes (Symonds and Moussalli, Reference Symonds and Moussalli2011). The models were then ranked by AICc and the one with the lowest AICc was selected as the best among all specified for each set. Then differences (Δi) between the AICc of the candidate models [AICc(min)] and the AICc of the best models [AICc(best)] were computed, and the Akaike's weight (wi) was calculated to obtain the relative likelihood of each model fitted of being the best among the whole set of candidates (Johnson and Omland, Reference Johnson and Omland2004).
In case of observing a high degree of uncertainty regarding the best model (when it is not strongly weighted; wi > 0·9), model averaging was computed on the subset of models whose cumulated weights (acc wi) were lower than 0·9, to estimate the magnitude of the effect of each predictor explaining most of the variation of the response variable. Moreover, to avoid redundancy, nested models of any model with a lower AICc value were excluded from the averaged subset (Grueber et al. Reference Grueber, Nakagawa, Laws and Jamieson2011). To assess the magnitude of the effect of each parameter in this subset, new Akaike weigths (w′i) were recomputed and parameter estimates, unconditional s.e. and unconditional confidence intervals (CIs) were calculated by averaging all models in the subset using full-model averaging method, weighted by the evidence of the respective models (models w′i) (Burnham and Anderson, Reference Burnham and Anderson2002; Symonds and Moussalli, Reference Symonds and Moussalli2011). In essence, this approach generates weighted average predictions, where the weights are the relative model probabilities. The continuous input variable (length) was standardized by subtracting its mean and divided by two times its standard deviation. Binary input variables (year, region of origin and sex) were coded using binary dummy variables and then were centred prior to multimodel inference analyses; this is strongly recommended when interactions are present because facilitates the interpretation of the relative strength of parameter estimates (Schielzeth, Reference Schielzeth2010; Grueber et al. Reference Grueber, Nakagawa, Laws and Jamieson2011).
All statistical analyses were conducted with the statistical software R version 3.0.1 (R Development Core Team; www.r-project.org/). The packages MuMIn (Bartón, Reference Bartón2009), arm (Gelman et al. Reference Gelman, Su, Yajima, Hill, Pittau, Kerman, Zheng and Dorie2015), AICmodav (Mazerolle, Reference Mazerolle2015) and the base package implemented in R were used.
RESULTS
Three species of myxozoans were recorded in M. hubbsi: K. rosenbuschi, M. meridionalis and a previously undescribed species of Fabespora (Table 2). Infections with A. merluccii or by M. baueri were not detected.
Geographical patterns of infection
Kudoa rosenbuschi was found in all samples, with higher values of prevalence in Patagonian samples; M. meridionalis was almost evenly distributed among all hake samples; while Fabespora sp. was only observed in Patagonian samples (Table 2).
Model validation plots (logistic regression quantile–quantile plots and partial residual plots) for every subglobal model in each set (see Tables 3 and 5) showed no major departures from model assumptions (graphs not shown).
All a priori selected biologically competing hypotheses (those considered as plausible explanations of parasite occurrences across samples in both regions), are identified and shown in Table 3 (13 for K. rosenbuschi and M. meridionalis, and five for Fabespora sp.).Within the set of models fitted for both K. rosenbuschi (K1–K13) and M. meridionalis (M1–M13), multiple models were well supported (within acc wi < 0·9, Table 3). Consequently, a model-averaging approach was used to incorporate information from these well-supported models in order to estimate parameter coefficients and standard errors (Table 4). In the case of Fabespora sp., for which region of origin was excluded from the analyses, the null model F2 lay within the subset of models whose acc wi < 0·9; therefore, the variables length and year, as well as their interaction, had no effect on the variability of occurrence of Fabespora sp. across hosts.
Parameter estimates (β̄̃), unconditional standard errors [ŝe (β̄̃)], 95% confidence intervals (CI) and relative importance (RI) of each predictor calculated under full model averaging methods of the occurrence of two myxozoan parasites (Kr – Kudoa rosenbuschi, Mm – Myxoproteus meridionalis) of Merluccius hubbsi in two regions of the South West Atlantic

L, total length; Y, year; R, region of origin; in bold predictor variable with significant effect on the occurrence of each myxozoan species.
Alternative candidate GLMs ranked according to AICc to evaluate the effects of predictor variables (L total length and PAT research fishing haul position) on the occurrence of three myxozoan parasites (Kr – Kudoa rosenbuschi, Mm – Myxoproteus meridionalis, Fsp. – Fabespora sp.) of Merluccius hubbsi in three sampling stations in Patagonian waters

ID, model identification; K, number of explanatory parameters in the model; AICc, second-order Akaike Information Criterion; Δi, AICc(best) – AICc(min); wi Akaike weight; acc wi, cumulative Akaike weight; LL, log-likelihood; over the dotted line models with acc wi < 0·9; in bold non nested models used for multimodel inference.
a Null model of each set.
b Subglobal model of each set.
Particularly for K. rosenbuschi, six (K1–K6) out of the 13 models fitted were well supported. Among them, only those providing non redundant information (K1, K2 and K4) were selected for multimodel inference based on its respective model w′i. These models consistently included region of origin as a predictor variable for explaining parasite occurrence and the inclusion of the other parameters explained no additional variation; in fact, they introduced additional uncertainty. Furthermore, model averaging showed that only region of origin had a CI 95% not including zero, indicating that this variable had a significant effect on the patterns of occurrence of K. rosenbuschi. On the other hand, for M. meridionalis, 8 (M1–M8) out of the 13 models were well supported, among which M1, M2, M4, M6 and M8 were not nested with models with lower AICc and, therefore, were used for multimodel inference based on its respective w′i (Table 3). After model averaging, year was the only predictor that explained a significant amount of variation of the occurrence of this parasite species (Table 4).
At a smaller spatial scale (within Patagonian samples caught in 2009), five biologically competing hypotheses were selected a priori as plausible explanations of parasite occurrences across samples (Table 5). Neither for K. rosenbuschi nor for M. meridionalis any predictor variable (year and length) explained the occurrence of these species as the null model (K′1 and M′1) of each respective set (K1–K5 for K. rosenbuschi and M1–M5 for M. meridionalis, see Table 5) lay within the subset of models with acc wi < 0·9. Notwithstanding, the results in Table 6 indicate that model F′1 with an Akaike weight of 0·71 is the best of five candidate models, indicating that fishing haul is the only predictor variable that explained the patterns of occurrence of Fabespora sp. within Patagonian samples caught in 2009.
Parameter estimates (β̄̃), unconditional standard errors [ŝe (β̄̃)], 95% confidence intervals (CI) and relative importance (RI) of each predictor calculated under full model averaging methods on the occurrence of Fabespora sp. on Merluccius hubbsi in three sampling station in Patagonian waters

Host sex had no effect on the occurrence of any myxozoan species within sample Patagonian 4 as the null model lay within the set of models whose acc wi < 0·9.
DISCUSSION
The available knowledge on the existence of recognizable stocks of M. hubbsi in Argentine waters, as well as on the occurrence of seasonal migrations between areas with different oceanographic characteristics provided the theoretical framework to test the suitability of myxozoans as biological markers of hake distribution. Their usefulness was corroborated, identifying them also as good metrics to detect geographical patterns at smaller geographic scales. Furthermore, with the finding of Fabespora sp., the list of myxozoans parasitizing M. hubbsi is extended to five species, increasing the potential of this group of parasites as BTs for hake stock discrimination.
Alatospora merluccii and M. baueri, previously recorded for M. hubbsi (MacKenzie and Longshaw, Reference MacKenzie and Longshaw1995; Kalavati et al. Reference Kalavati, Longshaw and MacKenzie1996), were not found in the present study. Indeed, both appear to have a southern distribution in the Atlantic Ocean, where they have been registered only around the Malvinas (Falkland) Islands (MacKenzie and Longshaw, Reference MacKenzie and Longshaw1995; Kalavati et al. Reference Kalavati, Longshaw and MacKenzie1996), suggesting that both myxozoans are coldwater species. These myxozoans constitute promising BTs for further studies on hake populations and migrations in that region.
Myxoproteus meridionalis was described by Evdokimova (Reference Evdokimova1977) from spores collected in the urinary system of M. hubbsi in the southwest Atlantic. Since then, it was recorded in this region in its type host (MacKenzie and Longshaw, Reference MacKenzie and Longshaw1995; Sardella and Timi, Reference Sardella and Timi1996, Reference Sardella and Timi2004) and in M. australis (MacKenzie and Longshaw, Reference MacKenzie and Longshaw1995). In the present study, this species showed no marked variations in prevalence among samples. MacKenzie and Longshaw (Reference MacKenzie and Longshaw1995) also found this species in hakes caught in the region of the Malvinas (Falkland) Islands, although at lower prevalence. Studies comprising more fish are necessary to evaluate whether this species could be considered a suitable BT to assess host population structure at a broader spatial scale, including hakes from the southern boundaries of its distribution. Meanwhile, M. meridionalis could not be considered a useful BT for discriminating stocks of M. hubbsi.
Kudoa rosenbuschi has been frequently reported from M. hubbsi (Sardella and Roldán, Reference Sardella and Roldán1989; MacKenzie and Longshaw, Reference MacKenzie and Longshaw1995; Sardella and Timi, Reference Sardella and Timi1996, Reference Sardella and Timi2004; Abollo et al. Reference Abollo, Novoa and Figueras2005; Whipps and Diggles, Reference Whipps and Diggles2006), showing an apparent geographical pattern. Indeed, MacKenzie and Longshaw (Reference MacKenzie and Longshaw1995) suggested a northern distribution of this species in samples from central to southern Patagonian waters, whereas Sardella and Timi (Reference Sardella and Timi1996) observed higher parasite burdens in southern areas when comparing hakes from Bonaerense and central Patagonian waters. Consequently, the central Patagonian waters seem to be the area where K. rosenbuschi is more prevalent. In agreement with previous findings, in the present study, K. rosenbuschi was found in all hake samples, but parasitizing a significantly larger proportion of fish from central Patagonian waters. San Matías gulf, the region where the higher prevalence has been recorded (Sardella and Timi, Reference Sardella and Timi1996), was not included in the present work. As proposed in earlier studies, K. rosenbuschi can be considered as a suitable BT for M. hubbsi. This parasite possesses the additional advantage that, due to its histozoic condition, with melanine pigmented pseudocysts easily detectable in skeletal muscle, its abundance can be estimated and used in quantitative comparisons. Nevertheless, it has been observed that K. rosenbuschi is more prevalent and abundant in juvenile hakes than in adults (Sardella and Roldán, Reference Sardella and Roldán1989). As an explanation of these patterns, Sardella and Roldán (Reference Sardella and Roldán1989) suggested that K. rosenbuschi could be an important cause of host mortality. However, as their sample was composed by hakes scattered in a broad region and caught during 3 years, this hypothesis remains unconfirmed. A decrease of prevalence with host age could be also explained by changes in the infection rate along host ontogeny, probably related to the development of immune responses, as has been suggested for K. allaria in Micromesistius australis (Agnew et al. Reference Agnew, Marlow, Lorenzen, Pompert, Wakeford and Tingley2003). The existence of ontogenetic variations in parasite burdens could impose limitations to the use of K. rosenbuschi as BT, since only hakes of similar size should be compared reliably. Nevertheless, in the present work, only hakes of similar size were included and, therefore, no effect of this variable was observed.
Fabespora sp. was only found in Patagonian samples; therefore it constitutes a promising BT for identifying Patagonian stock. Moreover, as samples were taken in the same year during two respective reproductive seasons of each hake stock, the differential presence of Fabespora sp. between regions (Bonaerense and Patagonian) suggests that these populations do not mix, at least not to any great extent.
Apart from large-scale differences, variations were also observed at a lower geographical scale. Indeed, among Patagonian samples (2009) Fabespora sp. infects a significantly higher proportion of fish in both southern samples, which in turn were mainly composed by female hakes (Patagonian 2 and 3) (Table 1). Skewed sex ratios are commonly observed for M. hubbsi during the reproductive season (summer) in the Patagonian area, a pattern attributable to local segregation produced by some reproductive behaviour that involves male-competition during fertilization, and to sequential movements of post-spawning females to deeper waters (Pájaro et al. Reference Pájaro, Macchi and Ehrlich2005; Macchi et al. Reference Macchi, Pájaro and Dato2007). Moreover, this pattern of sex segregation continues during feeding season (winter) forming ‘gender concentration nuclei’ along the Patagonian Shelf where M. hubbsi is also unevenly distributed by ages (Renzi, Reference Renzi1993). Consequently, it is difficult to discern whether the observed differences are attributable to host sex or to geographical origin.
Sex differences in parasite burdens may be caused by physiological or ecological processes, or by a combination of them (Reimchen and Nosil, Reference Reimchen and Nosil2001). On the other hand, the observed heterogeneity in parasite prevalence between hakes of Patagonian samples and, indirectly between host sexes may be the result of sex segregation during both trophic and reproductive seasons, with differential access to the endemic area of this myxozoan. However, Sardella and Timi (Reference Sardella and Timi1996) found no significant differences of the parasite burdens between sexes of hakes for none of the parasite species recorded. Similarly, in the present study, no differences between hake sexes were observed for any myxozoan species within sample Patagonian 4. As a consequence, the locality effect arises as the main explanatory variable of the observed differences between samples. Therefore, Fabespora sp. is proposed as a potential BT that would allow studying hake ecological features at various geographical scales (stock identity and intrastock temporal and spatial movements) in this region.
It is noteworthy that, at least for the promising parasite tags for hake population studies, namely K. rosenbuschi and Fabespora sp., no temporal patterns were observed for parasite occurrence. Indeed the year of capture, even after a period of 5 years, was not identified as a predictor of parasite burdens, indicating stability of parasite infections over long-term temporal scales. This temporal stability is a desirable feature for a BT and ensures repeatability.
Although the potential value of A. merluccii and M. baueri deserves further study, the present work confirms that the diversity of myxozoans harboured by hakes offers a set of suitable markers, which can provide evidence on stock discreteness at different spatial scales in the region. These species can be selected as BTs, either individually or in any combination, depending on the geographical extent of the study, constituting tools adaptable to the objectives of further research on hake population structure.
ACKNOWLEDGEMENTS
Thanks are extended to Marisa Lohfeldt, Claudio C. Burati (INIDEP) and Claudia V. Dato (INIDEP) for gently providing fish samples and to S.A. Barbini (Lab. Ictiología, IIMyC) for his valuable help with the statistical analyses. Special thanks are due to the 2 anonymous referees for their valuable comments that helped to improve this article.
FINANCIAL SUPPORT
Financial support was provided by grants from CONICET (PIP # 112-201101-00036) and ANPCYT (PICT 2012 # 02094) to DMPC, MMI and JTT.





