Post-genomic approaches to understanding the mechanisms of environmentally induced phenotypic plasticity

Living organisms rarely live in unchanging environments, and diurnal or seasonal fluctuations of abiotic environmental factors such as temperature,oxygen availability and salinity often dominate the life history strategies,physiologies and behaviour of organisms. However, it would be a mistake to regard exposed organisms as being entirely passive to the vicissitudes of life, since they frequently display a substantial capability to mitigate the direct effects of environmental change by physiological adjustments. These adaptive responses are particularly evident in species routinely exposed to challenge, and may be manifest in two ways(Cossins and Bowler, 1987), by maintaining normal levels of activity or of homeostatic potential (`capacity'adaptation), and by enhancing resistance to the potentially debilitating or lethal effects of environmental extremes (`resistance' adaptation). Exposure to environmental fluctuations during developmental stages may have life-long consequences, but adjustments made during juvenile or adult life are generally reversible. Both types of response can be regarded as prime examples of phenotypic plasticity.

Teleost fish from temperate, mid-latitudes have evolved within a strongly seasonal context, and because they are poikilotherms their tissues experience the full effects of the ensuing fluctuations in environmental temperature. Given that they have evolved in these varying environments they display some of the most powerful responses to diurnal and seasonal fluctuations of temperature. As a result cyprinid species such as the common carp, Cyprinus carpio, and centrarchid species such as the green sunfish of the mid-West of North America, have become favoured subjects for analysis of the mechanisms underlying environmentally induced physiological plasticity. The common carp is able to adjust to seasonal variations in temperature from<4°C up to >38°C, and water oxygen saturations down to just a few per cent of saturation. For changes in environmental temperature over the mid-range, these fish appear to adopt a compensatory strategy that sustains levels of metabolism and performance more or less constant despite the acute rate effects of temperature (Cossins and Bowler, 1987). This requires increases in the activity of enzymes involved in all manner of biological processes in animals exposed to cold and vice versa, brought about through changes in the cellular concentration and types of protein that are expressed(Hochachka and Somero, 2002). For adjustments outside the central range of temperatures, there may be other non-compensatory responses that can be regarded as leading to protection from the debilitating effects of extreme heat and cold(Cossins and Bowler, 1987).

Central to the whole body response to seasonal cold is the integration of adaptations in multiple, possibly all tissues, each of which may be manifest in many different ways. For example, in the common carp cold responses have been recorded in intestinal absorption through morphological and physiological changes to the absorptive epithelium (Lee and Cossins, 1988; Lee et al.,1991), to the performance of swimming musculature, through changes to the expression of Ca²⁺ regulatory and contractile proteins(Watabe et al., 1995) and to light-evoked properties of retinal horizontal cells(Cunningham and Hyde, 1995). Despite the fact that all parts of the body display thermal responses, our understanding of the underpinning mechanisms is generally fragmented, with only a few known genes that have been invoked as being involved in specific differentiated functions in each tissue. Also there is also no understanding of the common or even tissue-specific regulatory elements of responses that occur in tissues. Finally, there is little knowledge of the regulatory control that initiates and executes the process of adaptive change. Undoubtedly these responses require the coordinated activity of numerous genes and their encoded products, since temperature fluctuations affect all cellular and molecular processes, and there is a need for techniques to address not only the range of genes involved but their coordinated regulation. Understanding which genes are involved, in which tissues they are expressed, which levels of stress cause their induction and over what time course they occur, are prerequisites for advancing our understanding of environmental plasticity.

Whilst the analysis of candidate genes has provided some important information on particular aspects of the underpinning mechanisms, the responses have never been subjected to open screens for responding elements. Technological advances in recent years have provided powerful new techniques with which to undertake unbiased screens using transcriptome, proteome or metabolome data. In contrast to hypothesis-led analysis of candidate genes,these new screening methods potentially provide a system-wide assessment of response, to generate broad overviews of responses within which the role of defined biological pathways or processes and their underlying regulation can be interpreted. They also provide a valuable means of identifying new candidate genes that have not previously been identified through conventional hypothetico-deductive reasoning (Liang et al., 2004).

The full-blooded application of these techniques requires access to complete or near complete lists of genes, proteins or metabolites. Until recently this has limited attention to a small range of species for which complete genome sequences and gene lists are available, notably including fugu, mouse and human. However, recent work(Gracey et al., 2001) has demonstrated that the lack of sequence information is not an absolute barrier to progress; they generated 1600 cDNA probes for a microarray-based assessment of transcriptional responses of hypoxia responses in an estuarine gobiform fish species that is routinely exposed to environmental hypoxia. As a result,these contemporary techniques are now being applied to a much wider range of species, and to a range of environmental treatments(Cossins and Crawford, 2005)where they efficiently screen for genes displaying regulated expression. Studies of non-model species do not provide comprehensive genome coverage in that not all genes are sampled, but they tend to be genome wide, incorporating probes for all kinds of genes selected more-or-less randomly from the cDNA collection. Of course by focusing on genes displaying regulated changes in transcript expression, these methods do not address responses mediated at other levels, particularly including at the levels of protein turnover,protein post-translational modification, or epigenetic modification through histone manipulations.

To illustrate the way in which contemporary post-genomic techniques can advance our understanding of problems in comparative and integrative physiology, we describe here an analysis of the common carp, Cyprinus carpio, subjected to chronic cold exposure. This comprises the most extensive transcript screen to date of a non-sequenced vertebrate species responding to environmental challenge(Gracey et al., 2004). Here we focus on the responses of two of the seven tissues studied, namely intestine and muscle. We show how complex patterns of transcript responses can be interpreted within the framework of known physiology to provide new hypotheses, not only of the underlying physiology of phenotypic plasticity but also of the regulatory processes that control their activation.

The key requirement for pursuing a transcript screening approach is access to gene probes, either in the form of a collection of cloned cDNAs subjected to PCR amplification, or as oligonucleotides designed from the previously established gene sequence. In the absence of the latter we have selected∼14 000 individual clones from 14 directionally cloned, full-length cDNA libraries prepared from seven different tissues. These were subjected to single pass sequencing from the 5′ end. The resulting sequences were assembled into a minimally redundant set of contigs, and their identities were assessed by BLAST alignment against sequence databases. This resulted in a searchable database, carpBASE, which is available at<http://legr.liv.ac.uk>.

We amplified each clone by PCR and spotted the resulting amplicons onto glass slides using a robotic printer. These cDNA microarrays were then hybridised overnight with fluorescence-labelled cDNA prepared from RNA samples taken from the tissues of treated animals, and also at the same time with a reference RNA sample, common to all arrays, but labelled with a different fluor. Finally, hybridised arrays were scanned in both fluor channels and the two array images were quantitatively interpreted using proprietary software,which provides a relative measure of the binding of each cDNA sample to each spot on the array. Full details of these procedures can be found elsewhere(Gracey et al., 2004).

We used this basic microarray hybridisation procedure to monitor the gene expression changes that occur in the tissues of carp exposed to a regime of decreasing temperature, resulting in fish being held at 23, 17 and 10°C for up to 21 days. To control for natural day-to-day variation in gene expression we also analysed samples taken from 15 control specimens held at 30°C throughout and sampled on 3 different days over the full experimental period. We employed extensive replication throughout the experiment, and five or more individuals were sampled per time/temperature sampling point. We were particularly interested in how the response of different tissues contributes to the overall process of cold acclimation and so seven different tissues were dissected from each individual animal, resulting in ∼630 tissue samples. Following dissection the tissue samples were flash frozen and total RNA was extracted at a later date from a section of excised tissue cut from each frozen sample.

For this experiment we performed hybridisations of tissue RNAs against∼450 carp microarrays, each yielding a measurement of the relative mRNA level of each gene spotted on the array in each tissue sample versusthe level in a common reference RNA. Each RNA sample was hybridised against two arrays, each of which was labelled with a different combination of the two fluorescent dyes. From this large number of arrays, 374 passed a stringent quality control test and normalisation procedure that identified RNA samples that yielded congruent data on two arrays(Fang et al., 2003). The normalisation procedure removed spatial and intensity differences within and between arrays and provided an averaged measurement of the relative expression of each arrayed gene in each treatment RNA sample versus the common reference sample that was hybridised to each array. The analysis of expression across tissues thus comprised 374 individual values for each array probe,representing two repeated (dye-swap) analyses from each of 187 tissue RNA samples. The contribution of the reference sample to the expression data was removed by applying a `zero-transformation' procedure to all the data. This step simply involved dividing the expression level of a gene in each treatment sample by the average expression of that gene across the 15 control warm-acclimated animals for each tissue, thus providing an estimate of the relative expression of each gene in each cooled tissue sample versusits expression in the warm-acclimated control condition. The sheer amount of data in this experiment and the large amount of replication provided a high level of statistical precision (Gracey et al., 2004). Genes that exhibited a significant change in expression with cooling were identified using a popular signal-to-noise statistic (Tusher et al.,2001) and genes sharing similar patterns of expression were identified using a variety of pattern matching clustering algorithms(GeneSpring, Agilent, USA).

This procedure identified 3201 cDNA probes displaying significant changes in transcript expression at one or more time point during cooling relative to the warm-control specimens. We averaged the expression data for the arrayed cDNA probes that possessed identical BLAST identities since these likely represent redundant spots of the same gene, and explored the patterns in the whole dataset using principal component analysis. This revealed a common response to cold across all tissues, and comprised 260 unique genes that exhibited a significant and coherent change in gene expression in all of the seven tissues examined. Of these 260 genes, 221 shared homology with cDNA sequences that were already described in public databases. Of these 221 identifiable genes, all but eight were upregulated by cold, which is consistent with a general strategy of increasing the expression of proteins through transcriptional regulation as a means of compensating for the rate effects of cold, the increased protein expression offering additional enzymatic capacity to compensate for the rate-depressing effects of cold(Hochachka and Somero, 2002). They included genes that were associated with cellular homeostasis, notably including nucleic acid processing, cellular transport, protein catabolism,stress proteins and chaperones, metabolism, cell signalling and cell structure. Table 1 lists the genes with the largest fold-change values across tissues, most of which represent interesting candidate genes for further analysis. The list includes the Δ9-acyl CoA desaturase, a key gene that has previously been associated with cold responses and the increase in proportion of unsaturated fatty acids in the cold (Tiku et al.,1996; Polley et al.,2002). It also includes a glycine-rich RNA-binding protein similar to a cold-inducible RNA binding protein observed in the liver of Xenopus and other species(Nishiyama et al., 1997). Significantly, we show that it is induced in all tissues examined and we have also demonstrated that it is regulated during the hibernation cycle of the golden mantled ground squirrel (Williams et al., 2005). Brown and colleagues(Gerber et al., 2004) have recently described experiments with several members of the puf gene family in yeast, in which each protein binds a distinctive group of transcripts, probably via recognition sites in the 3′untranslated region of the transcripts. The carp protein identified in our work might also bind transcripts that are important in mediating the cold-stress response. Thirdly, the common response includes a gene with close similarity to uncoupling protein 3 (UCP3), which in mammals is expressed in skeletal muscle and brown adipose tissue, is upregulated by cold and is thought to participate in thermogenesis(Larkin et al., 1997; Ricquier and Bouilliard,2000). We find a wider tissue distribution and a substantial cold inducibility for this gene in all tissues examined(Gracey et al., 2004), but at present the physiological role of this gene is unclear apart from possibly being involved in energy metabolism, protection from reactive oxygen species damage and mitochondrial transport(Schrauwen and Hesselink,2002). A recent study(Jastroch et al., 2005) claims that UCP3 is restricted to skeletal muscle in the zebrafish, and the same is true in mammals (Schrauwen and Hesselink,2002), the difference with our analysis perhaps being related to different methods of transcript determination. This, together with the increase in expression on transfer to cold of both the ATP synthase and the ATP/ADP translocase, is consistent with the maintenance of ATP production despite the rate-depressing effects of cold exposure(Meerson, 1975).

The common response group includes many of the genes expected to be involved in basic functions of all cells, such as ATP synthesis, protein turnover, etc. Conversely they are distinct from genes whose expression underpins the differentiated functions of specific tissues. Interestingly,when we searched for orthologs of the carp common response genes in the yeast S. cereviseae, we discovered that the yeast orthologs for 25% of the carp genes were also found to be responsive to cold in a large yeast gene expression dataset (Gasch et al.,2000), indicating that a common set of genes is regulated by cold in both organisms This overlap suggests that the concept of a core response to environmental perturbation is generally applicable across a broad phylogenetic range, and that this group of genes has been conserved during long periods of evolution; as a result these genes may be used as diagnostic markers of cold stress in diverse tissues and organisms. Identifying a near-complete set of these conserved genes and understanding their distribution among species would offer powerful insights into how evolution has generated improved stress tolerance.

The remaining 1701 genes of known identity were differentially expressed in a single or as many as six tissues but not in all tissues. To address the enormous complexity of the tissue-specific responses, we clustered the expression data into 23 groups, each composed of 50–200 genes and with distinctive tissue-specific pattern of cold-regulated expression (see Fig. 1). Thus, cluster 2 represents genes whose expression was upregulated in 6 out of seven tissues,whilst cluster 5 included genes that were upregulated in just intestine only,and so on. Whilst these groups of genes showed similar tissue-specific patterns of regulation they comprise diverse genes that may participate in a variety of different processes. Identifying which of these cold-responsive genes is important can be achieved by profiling the cluster according the kinds of biological process involved(Hvidsten et al., 2003). For this we have used the conveniences of a controlled nomenclature developed for model species, namely the Gene Ontology (GO) database(Ashburner et al., 2000; Ashburner and Lewis, 2002). GO assigns genes, based on their established identity, into descriptive categories within three different domains: biological process, cellular function and molecular process. We have focused on the biological process domain and identified 24 categories that broadly represent all of the cold-regulated genes. We then developed an unbiased means of interpreting the genes within each cluster to identify the most prominent features by calculating the statistical probability that the representation of genes within a particular GO category was greater or less than expected from the number of genes interrogated on the array. This can be presented visually as a matrix of probabilities, called the GO-Matrix(Gracey et al., 2004).

Here we focus on the gene lists identified for the intestine to illustrate the way in which the GO-Matrix has been used to provide an unbiased interpretation. Previous work indicated that the intestinal mucosa of the carp is substantially remodelled in the cold by means of a doubling of intestinal wet mass, of rugal height and of mucosal surface area(Lee and Cossins, 1988). It also included the substantial upregulation of nutrient uptake across the mucosal border (Lee et al.,1991) and of ATPase activity of the basolateral membranes, as well as differential homeoviscous responses of mucosal and brush border membranes(Lee and Cossins, 1990). We previously suggested that mucosal growth in the cold involves a rebalancing between rates of enterocyte proliferation and losses through apoptosis, and this implies changes to underpinning mechanisms of either or both processes.

Of all the 387 responding genes identified as changing expression in the intestine, 88% were upregulated. Of these we found 172 possessing homology against genes in the sequence databases and with a GO annotation. For these we have used the main features of the GO-Matrix to guide the biological interpretation of the large-scale transcriptomics data. To illustrate this we focus on cluster 5, which largely consists of genes showing cold-induced upregulation in the intestine. The GO-Matrix analysis indicates that the following categories were significantly over-represented in the cluster:`Lipid metabolism', `Oxygen metabolism', `Transport', `Cell–cell signalling', `Cell adhesion', `Cell communication'. Only one GO category was under-represented, namely `Biosynthesis'. The predominant GO categories can then be interpreted within the framework of the known tissue-specific responses. Thus, the transport GO group(Table 2) includes 39 individual identified genes, some related to passive solute transport(Na/glucose transporter, Na/dicarboxylate transporter, aquaporin 9,facilitated glucose transporter, plasmolipin, neurotransmitter transport) and active transport (Na-K-ATPase subunit a-1). Others are related to mitochondrial transport of solutes, protons and electrons (cytochrome b, cytochrome c oxidase subunits VB, VIb and VIC-2,mitochondrial import receptor subunit, ATP synthase α, β andγ chains, phosphate carrier protein, uncoupling protein 3), to intracellular protein transport (ARF-related protein, ADP-ribosylation factor 1, GTP binding nuclear protein) or to lipid transport (Apolipoprotein A-1,A-1V, B-100, and E precursors, ATP synthase lipid binding mitochondrial protein).

Thus, the adjustments to the transport function of carp enterocytes during cold exposure can now be explored in the knowledge of these specific gene responses. Unfortunately not all gene probes identified through homology searching procedures were classified by the GO nomenclature and these were therefore not included in the GO-Matrix analysis. For the intestine, this comprised 84 genes with significant homology with gene databases but lacking any GO annotation, and also 141 genes lacking any identity. However, the list of transport-related genes greatly increases our knowledge of cellular processes undergoing regulation in cold-treated carp and this both supports and extends previous knowledge of mucosal and basolateral transporter systems.

An important use of microarray data is to provide some insights into the kinds of regulatory mechanisms that might initiate and direct the responses of transport and other genes. Table 3 includes all responding genes within the intestine that possessed GO annotation related to intracellular signalling, signal transduction and transcriptional regulation. This includes some high mobility group proteins, cell regulatory genes including a member of the DnaJ family, a series of kinases and phosphatases (MAPK-activated PK, PKC) and a series of transcription factors and regulators of transcription. Some of these genes are annotated as being involved in the regulation of cell cycle and these,together with several apoptotic genes, may lead to changes in the balance of cell proliferation and apoptosis. DnaJ is a large gene family with diverse functions in other cell types including suppression of cell death(Kurisu et al., 2003) and intracellular growth (Ohnishi et al.,2004). Similar arguments can be advanced for other genes known to be involved in intracellular regulation of proliferation and apoptosis,indicating a rich source of leads for exploring how the intestinal mucosa growth seen in the cold is initiated and controlled.

The breadth of microarray data and the ability to contrast responses between tissues or types of treatment enables the formulation of new hypotheses of the mechanisms of environmental response. For example, one of the most striking features of the carp cooling dataset(Gracey et al., 2004) was the discovery of a large cluster of genes that was repressed in white skeletal muscle tissue but unaffected in cardiac muscle(Fig. 2A). This shows that downregulation of this group of genes was specific to white skeletal muscle and is not common to other muscle types. GO profiling of this group indicated that it was significantly enriched for genes involved in cell motility and closer inspection revealed that the group included many genes associated with muscle contraction. These comprised genes encoding structural components of muscle fibers, for example, myosin light and heavy chains, as well as other genes, such as parvalbumin, with roles in handling sarcoplasmic Ca²⁺ levels and hence in muscle contraction and relaxation. Therefore, one interpretation of this expression signature is that it is indicative of an overall downregulation of the muscle contractile apparatus in response to cold.

In order to explore the mechanisms underpinning this response we searched for other sets of muscle genes that showed the opposite response in the cold,that is genes whose mRNA levels increased in skeletal muscle during cooling but which were unaffected in cardiac muscle. Fig. 2C shows one such cluster that met this criterion. Initial inspection failed to reveal a coherent biological theme regarding the function of cluster and therefore we inspected the genes one-by-one and searched the literature to identify those whose function was associated with muscle. This revealed two genes with a shared function in skeletal muscle, namely the muscle-specific ring finger protein 28(TRIM63) and the F box protein 32 (FBXO32)(Fig. 2B). These genes were first identified as essential components in the process of muscle atrophy in mouse (Bodine et al., 2001) and both are ubiquitin ligases that promote the ubiquitination and degradation of muscle proteins. This study showed that these genes were transcriptionally upregulated in response to reduced muscle activity and that ablation of either gene in mice promoted resistance to atrophy. Therefore, the increased expression of these genes in skeletal muscle of cooled carp together with the observation that many of their target muscle structural proteins are transcriptionally repressed is persuasive evidence that muscle atrophy may be activated at increasingly colder temperatures. The function of atrophy in the cold might be to reduce muscle mass through a process of reductive remodelling to balance for the reduced swimming activity observed in carp at colder temperatures. Note that cooling from 30°C to 23°C did not affect TRIM63 and FBXO32 transcript levels, though transitions to 17°C and then to 10°C induced substantial, but very transient,increases in transcript levels. This suggests that the response is induced only following more extreme cooling events, and that the response is graded according to the magnitude of the cooling imposed.

Previous work on muscle of adult cold-acclimated carp has emphasised the maintained performance of muscle in the cold, mainly through temperature-specific changes in the expression of myosin heavy chain isoforms(Tao et al., 2004; Watabe, 2002) and also for light meromyosin chains (Watabe et al.,1995). These cold-specific isoforms have properties that increase performance in the cold, perhaps due to possession of a more flexible protein structure. Cold also causes significant changes to the expression of myogenic regulatory factors during embryogenesis and development(Cole et al., 2004) and following the restoration of warm conditions in cold-acclimated adult carp(Kobiyama et al., 2000). These might be important in specifying the different muscle phenotypes in the cold. Our interpretation suggests that turnover of myofibrillar proteins is a major issue, either as part of the more general remodelling process or possibly to mediate a reduction in gross muscle performance under extreme cold conditions. This hypothesis remains to be tested at the level of protein expression.

Genomics-led approaches are by their nature large scale, and this demands the adoption of a new range of data processing and interpretation skills by investigators. However, depending on the technical approach adopted these approaches are unlikely to screen all possible genes; thus our microarray contained only 6033 non-redundant genes and even current versions of the zebrafish oligoarray contains just 17 000 of perhaps 25 000 genes. Given uncertainties regarding the size of the carp genome and how many genes its duplicated genome may contain, we estimate that the carp microarray represents just 12–24% of the entire complement of genes in the carp genome. Thus,the lists of responding genes generated by our analyses are likely to be just a small part of the overall picture. Even so, these genes are drawn from and thus sample the entire genome, and the resulting exploration of the partial transcriptomes can be highly informative of the full range of biological processes under investigation.

A second limitation of transcript screening is that it relates to just one level of biological regulation. There is a clear and well recognised need to validate key observations both at the level of transcript using RT-PCR and particularly at the level of protein or protein activity. Recently, there has been considerable discussion over the reliability of the microarray technique,not only with regard to the detection technique used(Drobyshev et al., 2003) and methods of handling noise that arises at all stages in the technique, but also to the statistical methods employed to discriminate significant changes in expression from non-significant effects through replication analysis. Whilst spurious and misleading observations can be obtained due to the complexity of the analysis, there is plenty of evidence in the literature that useful and meaningful results can be obtained provided that care is taken in processing the arrays, imaging and quantifying the arrays, and in the statistical analysis of the results (Irizarry et al.,2005; Larkin et al.,2005).

The number of potential leads generated by these open screening,system-wide approaches present new problems of interpretation. The first is the adoption of appropriate quality control and statistical procedures to identify which genes display significant responses. Importantly, this involves the use of intense replication in tissue sampling, the adoption of a suitable error model to estimate significance, and the use of an appropriate correction for the problem of false positives caused by performing a very large number of individual statistical significance tests. The statistical properties of microarrays are now generally well understood(Wit and McClure, 2004; Yang and Speed, 2002), and there is a considerable literature on the design of microarray experiments, on the downstream statistical processing and also the analysis of expression patterns through clustering. Our experimental design for the cold carp experiment incorporates a time series of exposure at three different levels of stress intensity, all sampled with replication to generate a very large data set. This enabled use of a powerful statistical approach with an error model encompassing nearly 400 different arrays, and this has enabled significant effects to be discerned at much low levels of fold-change.

Another problem is how to generate unbiased interpretations of the resulting complex patterns of gene regulation despite the temptation to follow well-trodden paths using genes whose properties are well understood. Given the large number of genes, and the potential for false positives, the permutations of interpretation that can be generated are endless, and interpretations thus need to be constrained both by the design of the experiment and by the methods of pattern analysis (Hvidsten et al.,2003). A less satisfactory alternative is arbitrarily to include some genes but not others, based on the judgement and knowledge of the investigator. On the other hand, some specific gene responses stand out as being unusual and call for further analysis. Good examples from our work include the need to understand the physiological significance of the RNA-binding proteins that are substantially upregulated in all tissues by environmental stress in carp (Gracey et al., 2004) as well as other species(Williams et al., 2005), or of the surprising non-muscle expression of myoglobin in carp tissues exposed to hypoxic treatment (Fraser et al.,2006). These two examples demonstrate the huge potential for identifying unexpected gene responses.

Transcriptomics using DNA microarrays provides a particularly tractable route to system-wide screening, particularly if genomics resources are readily available. Even if they are not, then techniques exist to generate sufficient resources within a moderately priced research project and over a period of a few months. This opens up an expanded range of species for which high throughput `-omic' approaches are tractable, including species of particular interest to the comparative physiology community. However, whilst transcript screening offers deep insights into only one level of biological regulation,they do not substitute for exploration of new candidate genes at the level of proteome or metabolome. Thus transcript screens are not an end in themselves,but act as a signpost pointing to further testable hypotheses. However,testing the role of these genes in generating specific, environmentally adaptive phenotypes will require application of gene manipulation techniques,which might be most easily addressed using more genetically tractable model species.

We thank Weizhong Li and Gregor Govan for skilled help, and the Natural Environment Research Council (UK) for funding.