A tale of two genes: divergent evolutionary fate of haptoglobin and hemopexin in hemoglobinless Antarctic icefishes

The isolated, frigid Southern Ocean has become well recognized as a vast natural laboratory where extreme polar conditions have wrought extraordinary evolutionary transformations. Today, its fish fauna is predominated by the five families of Antarctic notothenioid fishes that evolved in situ, forming an adaptive radiation and a rare marine species flock (Eastman, 2005). The success of these fishes in freezing waters owes to the evolution of a novel life-saving trait, the antifreeze glycoprotein (Chen et al., 1997), which inhibits the growth of ice crystals in blood and body fluids, thereby preventing organismal freezing (DeVries and Cheng, 2005). In contrast to such an adaptive fitness gain, the evolutionary transformations also included a no less impactful trait loss: the disappearance of the universal vertebrate oxygen transport protein hemoglobin in one of the families – the icefishes (Channichthyidae) (Sidell and O'Brien, 2006).

Icefishes are the only known vertebrates to be devoid of erythrocytes and hemoglobin (Hb) as adults. Instead, oxygen is transported as dissolved gas in the circulating plasma, at about only 10% of the concentration in the red-blooded notothenioids (Sidell and O'Brien, 2006). All 16 species of the icefish family lack Hb owing to the ancestral deletion of the β-globin genes in the Hb coding αβ-globin gene cluster (Cocca, et al., 1995), which occurred an estimated 8 to 9 million years ago (Near et al., 2006). In addition, six icefish species have further lost the expression of myoglobin (Mb), the key hemoprotein in oxidative muscles, in the heart. In contrast to the single origin of Hb loss, Mb loss was erratic, resulting from distinct mutational events in four separate occasions (Sidell and O'Brien, 2006).

Past research on the icefishes has greatly advanced our understanding of the mechanisms behind Hb and Mb loss, and the compensatory cellular, anatomical and physiological changes that followed (Sidell and O'Brien, 2006). However, little is known about the evolutionary impact of the hemoprotein losses on the network of protein partners that once supported their function. In the homeostatic balance of erythropoiesis [new red blood cell (RBC) formation] and eryptosis (destruction of senescent RBCs), partnering proteins are involved in the essential processes of recovery and disposal of free Hb and heme iron, which are otherwise cytotoxic because of their strong oxidizing properties (Chiabrando et al., 2011). As Hb and the RBCs that carry it headed towards extinction in the icefish ancestor, co-evolutionary changes very likely would have propagated across the supporting proteins. The fate of these supporting and partnering proteins has largely remained an unexplored part of the unique evolutionary history of the icefishes. Permanent absence of the hemoproteins would be expected to result in relaxation of selection pressure on the maintenance of their auxiliary systems (Lahti et al., 2009) and, as such, Hb loss (with Mb loss in addition) could be likened to the proverbial tip of the iceberg, with a much larger base of changes in partnering or supporting protein-coding genes.

Heme, with its iron center, is a strong oxidant, thus Hb will cause oxidative damage to surrounding molecules and tissues should it become free from its cellular carrier. The recovery and disposal of free plasma Hb and heme from RBC lysis owing to cellular senescence, infections or wounds depend on the direct action of two key plasma proteins, haptoglobin (Hp) and hemopexin (Hx) (Wicher and Fries, 2006). These bind with high affinity to their ligands, Hb and heme, respectively (Wicher and Fries, 2010). Therefore, they are excellent models for investigating the impact of relaxed selection from the Hb-less condition on the evolutionary fate of Hb-supporting proteins in Antarctic icefishes. Hp is well known as a scavenger of free Hb in vertebrates. In mice and humans, Hp has also been found to play various roles in the immune or inflammatory response. It modulates or regulates the balance of immune cells in mice (Arredouani et al., 2003; Huntoon et al., 2008) and stimulates angiogenesis in vasculitis (inflammation of blood vessels) in humans (Cid et al., 1993). Whether Hp in teleost fishes plays similar immune roles besides binding free Hb is currently unknown. Hx binds to free heme, which has sources beyond the Hb of RBCs. In the absence of RBC as the major source of free heme in icefishes, free heme could still be produced from other heme-containing proteins including mitochondrial cytochromes and microsomal cytochrome P450s, and other cellular hemoproteins including Mb in icefish species that express it. We hypothesize that the initial loss of Hb in the icefish lineage would commence a relaxation of selective maintenance of Hp and Hx functions, the extent of which would reflect the differences between the reservoirs of client proteins of these two scavengers. Here, we investigate the evolutionary fate of these two primary Hb and heme scavenger proteins in the Antarctic icefishes to gain new insight into the nature and extent of the impact of Hb loss on partnering molecules.

Most tissues for DNA and RNA extraction were taken from our existing inventory of preserved or frozen Antarctic fish samples. Fresh blood plasma for isolation of Hb-binding proteins was obtained for seven icefish species that we were able to catch in the West Antarctic Peninsula (WAP) waters during our field season in 2014. These included Chaenocephalus aceratus, Chaenodraco wilsoni, Champsocephalus gunnari, Chionobathyscus dewitti, Chionodraco rastrospinosus, Cryodraco antarcticus and Pseudochaenichthys georgianus. Blood plasma was also collected from two red-blooded notothenioids, Notothenia coriiceps and Akarotaxis nudiceps, to serve as positive controls for native Hp detection. Live fish were anesthetized with tricaine methanesulfonate (MS-222) (Western Chemical) at 1 g per 15 liters seawater, bled from the caudal vein using heparin-coated syringes, after which the fish were returned to the seawater holding tank to recover. The blood plasma was isolated after pelleting blood cells by centrifugation at 7000 rpm (4650 g) at 4°C for 20 min. Tissues were also dissected from several individuals of these species and either frozen in liquid nitrogen or preserved in 90% molecular grade ethanol at −20°C until use. All fish handling complied with the University of Illinois at Urbana-Champaign IACUC approved protocol 12123.

Plasma Hp and Hx are synthesized in the liver (Chiabrando et al., 2011), thus we extracted total RNA from liver samples of 15 of the 16 known icefish species, missing only Channichthys rhinoceratus, for which we lacked liver samples. RNA was also isolated from four red-blooded notothenioids, Dissostichus mawsoni, Pagothenia borchgrevinki and Notothenia coriiceps (family Nototheniidae), as well as Gymnodraco acuticeps (family Bathydraconidae), to serve as positive controls, for both RT-PCR amplification of Hp and Hx cDNA as well as northern blot analysis. Approximately 50 mg of tissue was homogenized in 1 ml TRIzol (Invitrogen) using 0.5 mm zirconium oxide beads and a Bullet Blender Storm^R (Next Advance), and RNA extraction from the homogenate followed the manufacturer’s (TRIzol) instructions. The final RNA pellets were dissolved in 50–100 μl of 0.5× TE (5 mmol l⁻¹ Tris–HCl, 0.5 mmol l⁻¹ EDTA, pH 8.0). RNA concentrations were estimated using absorbance at 260 nm (A260) and RNA purity by A260/A280 ratios (Epoch Take3 Microplate Spectrophotometer, BioTek Instruments). RNA integrity was verified by visualizing 1 μg of each sample run on a 2 mol l⁻¹ formaldehyde/1% agarose denaturing gel.

RT-PCR amplification for Hp and Hx cDNA was performed as an initial means to evaluate the transcriptional status of the two genes and the integrity of the protein coding sequence. Approximately 10 μg of each RNA sample was first treated with 2 units of DNase I (New England Biolabs) at 37°C for 10 min to degrade potential contaminating genomic DNA, then purified using E.Z.N.A. HiBind Mini Columns (Omega Bio-Tek). The DNA-free total RNA was primed with oligodT and reverse transcribed to generate the first-strand cDNA using SuperScript II reverse transcriptase (Invitrogen). To amplify Hp and Hx cDNA, gene-specific primers (Table S1) were designed using the genomic and transcriptomic resources available at the time, namely, the reference transcriptomes prepared for P. borchgrevinki (Bilyk and Cheng, 2013), D. mawsoni (Chen et al., 2008) and the icefish C. aceratus (Shin et al., 2012).

For Hx, a pair of UTR primers, Hx_UTR_F and Hx_UTR_R successfully amplified Hx cDNA from both red-blooded notothenioids and all 15 icefishes. For Hp, the UTR primers Hp_UTR_F and Hp_UTR_R successfully amplified full-length Hp cDNA from four icefishes but not the others (Table S1). A priori, it was unknown whether the null result was due to absence of Hp transcription or primer mismatch. Pairing the Hp_UTR_R primer with a forward primer Hp_ATG_F that anneals to the start of the protein coding sequence (the first 27 nt inclusive of ATG start) produced an amplicon covering nearly full-length CDS (coding domain sequences) for three additional icefish species. This suggested that divergence in the 5′ UTR sequence in various icefishes may have contributed to the null amplification. We thus attempted 5′RACE to obtain species specific 5′ UTR sequences to design the UTR forward primer but without success, suggesting absence or extremely low level of transcription not amenable to amplification. Thus in total, we obtained full-length protein CDS from four icefish species, and nearly full-length CDS for three other icefishes. Amplified Hp and Hx cDNA were sequenced with BigDye Terminator v.3 Cycle Sequencing chemistry (Applied BioSystems) by direct sequencing or following cloning into the pGemT_easy vector (Promega).

Hp and Hx protein coding sequences were investigated for evidence of changed selective pressure upon Hb loss in icefishes. We used RELAX (Wertheim et al., 2014), a general hypothesis test using a codon-based phylogenetic framework to assess whether selective constraint became relaxed or intensified gene-wide in test species compared with reference species. The complete coding sequences of Hp, Hx and mitochondrial ND2 genes for a subset of teleost fishes were obtained from Ensembl to serve as the reference set. Additional notothenioid Hp and Hx sequences came from the Parachaenichthys charcoti genome (Ahn et al., 2017) and the Eleginops maclovinus reference transcriptome (Bilyk et al., 2018). We then performed six comparisons to detect molecular evolutionary signals: (i) all icefish branches against all red-blooded Antarctic notothenioid branches, (ii) all icefish branches against reference teleosts, (iii) all red-blooded Antarctic notothenioid branches against reference teleosts, (iv) the single ancestral branch to icefish against all other icefish branches, (v) the single ancestral branch to icefish against all red-blooded notothenioid branches, and (vi) the single ancestral branch to icefish against reference teleosts. Table S2 lists the species and accession numbers of sequences used in this analysis. Only the four icefish species with full-length Hp coding sequences were included in this analysis. To reconstruct phylogenetic relationships of these species needed for the RELAX analyses, we used the ND2 sequence rather than Hp or Hx, as the latter do not possess enough informative sites to resolve the relatedness of icefishes. We aligned the complete ND2 gene sequences with codon constraint using MUSCLE (Edgar, 2004), and evaluated for models of nucleotide substitution using jModelTest (Posada, 2008). The best fit model GTR+I+G was implemented in Bayesian phylogenetic analyses using MrBayes v.3.2.6 (Huelsenbeck and Ronquist, 2001). We ran the Markov chain Monte Carlo simulation for 4 million generations with four chains and sampled every 100 generations. We then used Tracer 1.4.1 (Drummond and Rambaut, 2007) to examine the trace files to ensure the chains reached convergence, and discarded the first 25% of trees as burn-in. The consensus tree (Fig. S2) was used as the phylogenetic framework for RELAX analysis. Translated protein sequences for Hp and Hx were aligned using MUSCLE, and the sequences were manually trimmed to remove gap regions and premature stops. The codon aligned nucleotide sequences were then analyzed using RELAX from the HyPhy software package (Pond et al., 2005) implemented on the Datamonkey Adaptive Evolution Server (http://datamonkey.org/; Weaver et al., 2018).

We utilized northern blot to provide a tissue-wide view of the transcriptional status of Hp and Hx in liver among the icefishes, and as a means to validate loss of transcription where RT-PCR amplification produced no amplicon. Northern blot analysis was carried out on the 15 icefishes with available total liver RNA, and four representative red-blooded notothenioids as positive controls. Eight micrograms of total RNA was resolved on a 1.1% agarose formaldehyde gel at 110 V for 65 min, and imaged. The gel was gently shaken at room-temperature sterile-distilled water for 30 min to remove excess formaldehyde, then vacuum-blotted (Amersham VacuGene, GE Health Sciences) onto Hybond N membrane (Amersham, GE Health Sciences) and UV crosslinked (Stratalinker, Agilent Technologies). The membrane was prehybridized at 55°C in PerfectHyb (Sigma-Aldrich) for several hours, then hybridized at the same temperature overnight to P³²-dATP (Perkin Elmer) labeled probes of Hp or Hx. Probes were synthesized using random heptamer priming of templates consisting of equal amounts of Hp or Hx cDNA derived from the red-blooded notothenioid D. mawsoni and the icefish C. aceratus. The hybridized blot was washed with increasing stringency, to 0.1× SSC (15 mmol l⁻¹ sodium chloride, 1.5 mmol l⁻¹ sodium citrate)/0.5% SDS (sodium dodecyl sulfate) at 55°C for Hp and 57°C for Hx. The washed blot was autoradiographed on a phosphor storage screen (Kodak) overnight and scanned on a STORM^R 860 phosphoimager (Molecular Dynamics) to obtain the image. Between the two gene probes, the blot was stripped of hybridized probes by immersion in 0.1% SDS that was heated to boiling.

To characterize the control region of the Hp gene, we sequenced select Hp-positive large-insert DNA genomic clones isolated from BAC (bacterial artificial chromosome vector pCC1BAC, Epicentre) libraries of the icefish C. aceratus (Detrich et al., 2010), the red-blooded Antarctic notothenioid D. mawsoni (Nicodemus-Johnson et al., 2011) and the basal temperate notothenioid Eleginops maclovinus for evolutionarily relevant comparison. The D. mawsoni and E. maclovinus BAC libraries had been constructed in our lab for use in various studies, whereas the C. aceratus library was available through Children's Hospital Oakland Research Institute (CHORI) BAC Resources Center. The macroarray filters of the BAC libraries were screened for the Hp loci using P³²dATP-labeled probes generated from the C. aceratus or D. mawsoni Hp cDNA using the same hybridization process described above but washed less stringently, to approximately 50°C. Putative positive BAC clones were isolated from the D. mawsoni and E. maclovinus libraries, and those of C. aceratus were purchased from CHORI. These were further verified by Southern blot probing of the NotI digested recombinant BAC plasmid DNA. The positive BAC clones were found to form a single contig group on fingerprinted contig (FPC) analysis (Soderlund et al., 2000), thus there is a single Hp locus represented by any of the BAC clones within the group. The plasmid DNA of one selected positive clone for each species was electroporated into the TransForMax EPI300 Escherichia coli strain (Epicentre), which was induced to produce the plasmid at high copy number to generate sufficient DNA for sequencing. Sequencing libraries were constructed using the Roche GS Titanium Nextera DNA Sample Prep Kit and sequenced on a Roche/454 GS FLX+ (Roche/454 Life Sciences) at the Roy J. Carver Biotechnology Center at the University of Illinois at Urbana-Champaign. The sequence reads were assembled for each species using Newbler v2.7. The sequence contig containing the Hp gene was analyzed for the presence and integrity of proximal promoter elements within approximately 1 kbp upstream of the start codon using JASPAR (Khan et al., 2017).

The sequenced Hp loci were used to design a 3′ flanking region primer (Hp_3Flank_R), which was then paired with a forward primer (Hp_E4_F) to amplify much of the large terminal exon (exon 4) and part of the downstream 3′ flanking sequence of the Hp gene from genomic DNA of all 16 icefish species and two red-blooded control species, D. mawsoni and G. acuticeps. This provided a means to verify the continued presence of the gene within the genomes of the icefish species where we could not amplify cDNA products. Genomic DNA was isolated from liver or muscle using standard Tris/EDTA/SDS tissue lysis, phenol/chloroform extraction and final ethanol precipitation. PCR amplified products were treated with SAP/ExoI (shrimp alkaline phosphatase/Exonuclease I) and then sequenced with BigDye Dye Terminator v.3.1 chemistry (Applied BioSystems) by direct sequencing or after cloning into pGemT_easy vector (Promega).

Fresh native plasma collected for seven icefish species were used for immune detection of physiological levels of circulating Hp protein. Anticipating low physiological Hp levels, we also enriched the Hb-binding protein fraction by affinity chromatography. Hb of the plentiful WAP species N. coriiceps was used as bait to isolate Hb-binding proteins. Spun N. coriiceps erythrocytes were washed twice with notothenioid phosphate buffered saline (PBS; 12 mmol l⁻¹ NaH₂PO₄, 86 mmol l⁻¹ Na₂HPO₄, with osmolality adjusted with NaCl to 550 mOsm, the physiological level of native notothenioid plasma), then lysed in five volumes of an ice cold, hypotonic solution of 10 mmol l⁻¹ Tris (pH 8), and centrifuged at 16,000 g for 10 min at 4°C to pellet cellular debris. The supernatant was chromatographed by gravity flow on a 1.5×20 cm column of a Sephadex G-75 (GE Healthcare) in 10 mmol l⁻¹ Tris (pH 8), and the visibly red Hb-containing elution fraction was collected. This purified Hb was coupled to CNBr-activated Sepharose™ 4B (GE Healthcare) following the manufacturer’s instructions, and then packed in a Polyp-Prep column (BioRad). To capture Hb-binding proteins inclusive of Hp, 1.5–2 ml plasma from each species was diluted 3-fold with notothenioid PBS, and run through the column three times by gravity flow. After thoroughly washing with 30 column volumes of notothenioid PBS, the bound proteins were eluted with 5 ml of 8 mol l⁻¹ urea. The eluate was then dialyzed overnight in 10 mmol l⁻¹ Tris (pH 8) at 4°C and the dialyzed proteins were lyophilized.

To detect the presence of Hp in native plasma, and in the affinity column enriched Hb-binding protein fraction, we used an anti-Hp polyclonal antibody in western blot analysis. The anti-Hp antibody was custom produced (by GenScript) in rabbit against several in silico determined antigenic epitopes from conserved regions of notothenioid Hp sequences we obtained for icefishes. Synthetic peptides of these epitopes were made with a terminal cysteine added for conjugation to keyhole limpet hemocyanin. Antibodies generated with the epitope DKVTPIPLPERGQDC provided the strongest binding and sensitivity in subsequent tests. To first test the specificity and sensitivity of the anti-Hp antibodies, we expressed and purified recombinant Hp of the icefish Chionodraco hamatus, which was also used as an icefish Hp positive control in western immunoblotting, as well as for determining appropriate antibody dilution factors. The C. hamatus Hp cDNA was cloned into the pET26b+ expression vector upstream of the HisTag and transformed into Rosetta2(DE3)pLysS host cell (Novagen, Millipore-Sigma). A verified recombinant cell clone was grown at 16°C in lysogeny broth containing 50 µg ml⁻¹ kanamycin and 34 µg ml⁻¹ chloramphenicol to an OD₆₀₀ of 0.6, at which point protein expression was induced by adding IPTG (isopropyl β-d-1-thiogalactopyranoside), followed by further growth of the culture for 4 h. Cells were harvested by centrifugation and lysed in BugBuster^R protein extraction reagent (Novagen, Millipore-Sigma). Hp was extracted from inclusion bodies and then purified using His-bind column chromatography (Millipore-Sigma) following the manufacturer’s instructions. The purified Hp was eluted in buffer containing 1 mol l⁻¹ imidazole and 6 mol l⁻¹ urea, and the concentration was measured using a Bradford protein assay kit (BioRad).

Western blot analyses of Hp were carried out on samples of native plasma samples and the lyophilized, affinity column enriched Hb-binding protein fraction of the two representative red-blooded species N. coriiceps and A. nudiceps, and the seven icefish species with available fresh plasma, and the icefish C. hamatus recombinant Hp as positive control. Protein samples were resolved on precast 8–12% gradient SDS-PAGE gels (Invitrogen), and electrophoretically transferred to Immun-Blot PVDF membrane (Bio-Rad). The membrane was then shaken gently in 5% gelatin (Sigma) for 1 h at room temperature to block non-specific sites, then incubated for 2 h in PBST (PBS with 1% Tween 20) with a 1:1000 dilution of the primary anti-Hp antibody (stock concentration 0.864 mg ml⁻¹), followed by washing three times with PBST. The washed blot was incubated for 1 h with the secondary goat anti-rabbit antibody diluted 1:20,000 in PBST from a 0.8 mg ml⁻¹ stock. After washing with PBST, the membrane was treated with enhanced chemiluminescence (ECL) western blotting substrate (Pierce), and the signal was imaged with a LI-COR C-DiGit Imaging System.

We examined the icefish Hx sequence features and modeled their protein structures to assess the integrity of structural components essential for Hx function. A mammalian (rabbit) Hx sequence (Morgan et al., 1993) was used as reference to identify the distinctive sequence features and functional residues of Hx in notothenioid sequences. Using the X-ray crystal structure of this Hx (Protein Data Bank accession number: PDB|1QHU|A chain, doi:10.2210/pdb1que/pdb) (Paoli, et al., 1999) as template, Hx protein structural models of two icefishes, C. wilsoni and P. macropterus, and the red-blooded N. coriiceps were obtained by homology modeling using the SWISS-MODEL server (Waterhouse, et al., 2018) to assess global and local folds.

We successfully amplified full-length or nearly full-length Hp cDNA by RT-PCR for seven of the 15 available icefish species and four red-blooded control species. The sequences have been deposited in GenBank under accession numbers MH548902–MH548912. Full-length cDNAs were obtained for four icefish species: Pagetopsis macropterus, C. aceratus, C. wilsoni and C. hamatus. The C. wilsoni and C. hamatus Hp cDNAs encode an intact protein of 315 amino acids, whereas C. aceratus and P. macropterus cDNAs contained independent nonsense mutations leading to premature termination at amino acid positions 24 and 133, respectively (Fig. 1). Potential 5′ UTR sequence divergence and/or extremely low or absence of Hp transcription precluded amplification of full-length cDNA in the other icefish species. The use of an inner primer that anneals to the first nine codons (27 nt) of the signal peptide to pair with the 3′ UTR primer succeeded in amplifying partial (nearly full length) CDS for three additional icefish species, Neopagetopsis ionah, Dacodraco hunteri and Chionodraco myersi (Fig. 1). These nearly full-length cDNAs serve as evidence that the Hp gene in these three species remained transcriptionally active. Interestingly, N. ionah Hp contained a 12 nt deletion (residues 169–172) without resulting in a reading frame shift, while D. hunteri Hp sustained the same position 24 nonsense mutation as C. aceratus, and additionally a downstream frame-shift mutation at positions 52–54 (Fig. 1, Fig. S1).

In contrast to Hp, we successfully RT-PCR amplified Hx cDNA from all available (15 of 16) icefish species and the four red-blooded controls, encoding a protein of 426 to 443 amino acids depending on species (Fig. 2). These Hx sequences have been deposited in GenBank under accession numbers MH546081–MH546099. Readily amplifiable icefish Hx cDNAs indicate sustained transcription of this gene throughout the icefish family. The Hx protein coding sequences show no apparent molecular defect (Fig. 2), indicating continued integrity of the Hx gene in the icefishes.

With the phylogenetic framework for the target notothenioids species and other teleost fishes (Fig. S2), the RELAX test identified a significant relaxation in selective pressure on the Hp sequence of both the Hb-less icefishes and the red-blooded notothenioids relative to the reference teleost fishes. In both comparisons, the estimated intensity of selection parameter, K, was found to be below 1, indicating that selection strength has been reduced along the test branches (Table 1). However, direct comparison between the icefishes and the red-blooded Antarctic notothenioids failed to find a significant change in selective pressure, as did all of the comparisons against the ancestral icefish branch. In contrast, Hx showed signs of intensification in selective pressure (K>1) in both the icefishes and the red-blooded Antarctic notothenioids when compared against the other teleosts (Table 2).

To ascertain whether technical issues (primer mismatch, non-optimal RT and PCR conditions, etc.) solely contributed to unsuccessful amplification of Hp cDNA in some of the icefish species, we performed a northern blot to broadly assess liver Hp transcription (Fig. 3). The results showed that among the 15 examined icefish species, Hp mRNA was weakly detected in five of the icefishes: P. macropterus, C. aceratus, C. wilsoni, C. hamatus and C. rastrospinosus (Fig. 3B), three of which had produced full-length Hp cDNA on RT-PCR amplification. The transcript abundances, however, are miniscule compared with those of the red-blooded notothenioids. Contrasting with Hp, Hx transcription appears qualitatively robust in the icefishes, albeit the levels based on hybridization intensity are highly variable, but the variability was similarly observed in the red-blooded species (Fig. 3C).

Sequencing results for the Hp genomic loci in E. maclovinus, D. mawsoni and the icefish C. aceratus are deposited in the NCBI SRA (accession numbers SRR5878014–SRR5878018), and the Hp gene structure is diagrammed for all three species in Fig. 4A. The complete gene was assembled into a single contig for E. maclovinus and D. mawsoni, but was split between two contigs in the icefish C. aceratus. However, this split did not preclude obtaining the complete protein coding sequence, or intron/exon boundaries. The Hp CDS from spliced exon sequences of the gene from the C. aceratus locus assembly corroborates the occurrence of the premature stop codon at amino acid position 24 found in the cDNA sequence (Fig. 1). We attempted to amplify the Hp gene from genomic DNA of the icefish species for which we could not obtain cDNA, but this was unsuccessful. The 5.7 kbp Hp genomic sequence of D. mawsoni is quite long, primarily owing to the large intron 2 (Fig. 4A), but we were able to PCR amplify this sequence from genomic DNA by long-distance PCR (results not shown). The failure to amplify the icefish Hp genomic sequences suggests that intron 2 might have further expanded, perhaps at the repeat-rich region that had led to the split of the gene in our assembly of the C. aceratus Hp genomic locus from BAC clone DNA. We therefore utilized the 3′ flanking sequence of the Hp locus of C. aceratus to develop a reverse primer and paired it with a forward primer at the start of exon 4, and successfully amplified a large part of the Hp gene for all 16 icefish species in the family. Readily amplifiable Hp genomic sequences indicated that the Hp gene (or at least the bulk of the gene amplified) persists in the genomes all 16 icefishes. The partial sequences are shown in Fig. 4B, and have been deposited in GenBank under accession numbers MH546100–MH546112. These partial Hp sequences revealed another distinct premature stop (position 39 of the partial sequence; Fig. 4B) in Pseudochaenichthys georgianus and Cryodraco antarcticus at amino acid position 174 in the full-length protein (Fig. 1). In addition, we discovered a mutated polyadenylation signal sequence, AATAGA instead of the canonical AATAAA, in five of the species (Fig. 4C).

Using the Hp genomic loci of E. maclovinus, D. mawsoni and the icefish C. aceratus, we searched for control element motifs within the proximal 1 kbp sequence ahead of translational start site and located a logically placed TATA box at approximately 80 bp upstream from the ATG start in all three species. Additionally, putative cis-promoter sites for transcriptional factors including HNF1, HNF4, CEBPB/D and DBP found to be present in the two red-blooded species are conserved in the icefish (Fig. 5).

Expecting native Hp levels would be low in icefishes, a recombinant C. hamatus Hp protein was generated to confirm the specificity and sensitivity of the Hp antibody, and the titration western blot showed that the antibody can detect Hp down to as low as 0.2 ng µl⁻¹ (Fig. S3A). Subsequent western blot analysis on native blood plasma (Fig. 6A,B) revealed a single strong immunopositive band in the control red-blooded notothenioids N. coriiceps and A. nudiceps at an estimated size of 37 kDa, close to the average molecular mass of approximately 35 kDa calculated from translated Hp sequences. In contrast, western blot of native blood plasma from the icefishes failed to detect Hp at their natural plasma concentrations.

The failure to detect Hp in icefish blood plasma could reflect an absence of expression or natural expression at levels below the 0.2 ng µl⁻¹ sensitivity of our anti-Hp antibodies. Hb from red-blooded N. coriiceps was used as bait in an affinity column for testing the Hb binding ability of Hp in the blood plasma of the Antarctic icefishes if present. This affinity column also serves to concentrate the icefish plasma Hp from physiological levels, where expressed, providing a highly sensitive assay for its presence. Using a Hb-binding protein fraction enriched by affinity chromatography, the presence of Hp was again confirmed in the red-blooded notothenioids (Fig. 6C,D), and also in one of the seven icefishes, C. wilsoni, showing a putative 40 kDa Hp protein (Fig. 6B,C, Fig. S3B). The enrichment by Hb affinity column was approximately 300-fold, which by definition means the C. wilsoni Hp is capable of binding Hb, but is expressed at extremely low levels.

We examined the canonical sequence features and functional residues, as well as protein structure models of translated Hx protein sequences of the notothenioid fishes to assess functionality. The well-studied mammalian Hx consists of two distinctive homologous Hx domains connected at 90 deg to each other by a short linker, with each domain comprising four sequence repeats that form a four-bladed β-propeller structure. Three pairs of cysteine residues in each domain form three disulfide bridges, and a pair of histidine residues, one in the linker and the other in repeat 1 of the C-terminal domain, coordinate the heme ligand and stabilize it in the binding pocket between the two domains (Paoli et al., 1999; Piccard et al., 2007; Shrimal and Gilmore, 2013). These sequence features and functional residues are conserved in both the red-blooded and white-blooded notothenioids as indicated in the sequence alignment (Fig. 2) with additional details for three representative species selected for Hx structural modeling (Fig. S4). Using the mammalian (rabbit) Hx X-ray crystal structure (1HQU chain A) (Paoli et al., 1999) as template, homology modeling of the translated Hx sequences of a representative red-blooded notothenioid N. coriiceps and the two icefishes C. wilsoni and P. macropterus showed high degree of conservation of global structure and the two four-bladed propeller domains (Fig. 7). The GMQE (global model quality estimation; range 0 to 1) scores of all three notothenioid models are 0.71, indicating high reliability of the models. QMEAN (qualitative model energy analysis) scores are −1.45, −1.37 and −1.70 for N. coriiceps, C. wilsoni and P. macropterus, respectively, indicating good agreement between the modeled structures and experimental structures as they are much higher than the −4.0 cut-off score for low-quality models (Benkert et al., 2011). These results are consistent with the notothenioid Hx cDNAs being translatable into a functional protein.

The evolution of Hp and Hx in vertebrates as scavengers of cytotoxic free Hb and heme was intimately tied to the appearance of Hb, emerging in fishes as Hb became present at high corpuscular concentrations (Wicher and Fries, 2010). The ancestor of the Antarctic icefishes suffered genetic loss of this important oxygen transport protein, yet survived and diversified into 16 species (Kock, 2005). This trait loss begs the question of how it in turn affected the fate of partnering protein genes that once supported erythropoiesis and eryoptosis, particularly the Hp and Hx proteins with dedicated Hb and heme binding function. This study provides the first examination of this evolutionary domino in the Antarctic icefish family, and the outcome from the release of selective pressure due to the absence of RBCs and Hb. The contrasting fates we observed for the two candidate genes offers insight into the impact of Hb loss that may fall on the network of its supporting genes.

The Hb-null state of the entire icefish clade means that the selective pressure to maintain functional Hp as an Hb-scavenger would have been removed at the origin of the icefishes, with the expectation that the trait would continue to degenerate and potentially become extinct. However, Hp exhibits a peculiar pattern of persistence and loss across vertebrate taxa. Chicken (bird) and western clawed frog (amphibian) lacked the Hp gene, with a different scavenging protein (PIT54) replacing Hp as the free Hb scavenger in chicken (Wicher and Fries, 2006). Hp is present and expressed at high levels in cartilaginous fishes, but its Hb binding affinity is low or nil, suggesting that Hp may serve other, as yet unknown, non-Hb scavenging function in these fishes (Redmond et al., 2018). In teleost fishes, the Hp gene is present in various species, and native serum of the Japanese pufferfish binds Hb (Wicher and Fries, 2006), but trout Hb does not, again suggesting it persists in trout for a non-Hb binding role(s) (Redmond et al., 2018). With only two teleost species in which Hb-binding has been assessed and producing contradictory results, it remains unclear whether Hp in general has other non-Hb binding functions in teleosts in general, or whether it is species specific, and what the function(s) may be. In this study, we showed that the Hp gene is present in the basal notothenioid E. maclovinus and across the derived Antarctic species, both red- and white-blooded. Further, functional results indicate an Hb binding by a plasma protein (very likely Hp) in the red-blooded Antarctic notothenioids N. coriiceps and A. nudiceps (Fig. 6B) and the icefish C. wilsoni (Fig. 6D). These establish that the Hp genotype and Hb-binding function is a shared plesiomorphic trait in the Notothenioidei suborder. If the only or predominant role for Hp is Hb binding in the notothenioids, then the loss of Hb at the origin of the icefishes would have removed its major, if not the sole, client protein, relaxing selective constraint for its maintenance.

The hypothesized loss in selective pressure on Hp maintenance was corroborated by the near silencing of Hp transcription and translation across the icefishes we examined. With the less sensitive but more global detection by northern blot, Hp transcription in icefishes has greatly declined in comparison with the red-blooded relatives and appears non-existent in several species (Fig. 3B). With repeated optimization, the more sensitive RT-PCR succeeded in amplifying Hp cDNAs from seven of the 15 available icefishes, indicating that Hp transcription exists at low levels in these seven species (Fig. 1). The inability to amplify Hp cDNA by RT-PCR or detect mRNA by northern blot for the other eight species indicatesthat either the transcription level is below the detection limit of these two methods, or that these species may have lost Hp transcription entirely. The mutated polyadenylation signal sequence in five species – C. rhinoceratus, C. aceratus, C. dewitti, C. atkinsoni and C. rastrospinosus (Fig. 4C) – likely contributes to mRNA instability, obviating ready detection.

Investigations of the translational status of the Hp gene provided an even greater contrast. Only one icefish species, C. wilsoni, showed a protein of expected size for Hp on western blot (Fig. 6C,D) after its Hb-binding protein fraction in the plasma had been greatly enriched (∼300-fold) for Hp using Hb-affinity column chromatography. Thus, physiological Hp protein levels must be extremely low even if expressed, and below the threshold of detection (0.2 ng μl⁻¹) of our antibodies (Fig. S3A). No Hp could be immuno-detected in the other six icefish species with available fresh plasma, whether in the native plasma (Fig. 6A,B) or affinity column enriched Hb-binding protein fraction from a large volume of plasma (Fig. 6C,D). This shows a paucity of functional protein product at physiological conditions either from translation levels too low to measure, or inability of any synthesized Hp to bind Hb, or no Hp protein synthesis at all.

Although terminal evidence of relaxed selective pressure, in the form of loss-of-function mutations, occurred only in the icefishes, tests for shifts in selective pressure detected signs of relaxation in both the icefish and the red-blooded Antarctic notothenioids (Table 1). Additionally, direct comparison against the red-blooded Antarctic notothenioids failed to find a significant signature of further relaxation specific to the icefishes. One possible reason for this shared relaxation in selective pressure is a potential lower incidence of free plasma Hb, even in the red-blooded species, owing to much lower hematocrits and mean corpuscular Hb concentrations in Antarctic notothenioids compared with temperate fishes (Wells et al., 1980), and with the expected slower turnover of RBCs given their extreme low body temperature. But it is only among the icefish species in which the need for Hb binding had been fully removed that nonsense mutations were capable of reaching fixation. Ultimately, if there is only relatively weak further relaxation in selective pressure with the origin of the icefishes than this may not be readily detectable using existing approaches.

In marked contrast to Hp, northern blot showed that Hx is clearly transcribed across all 15 of the tested icefishes, some at prominent levels (Fig. 3C). The Hx mRNA levels are variable among icefish species, but variations are also observed in the red-blooded notothenioids. Whether they resulted from differing levels of physiological demands for Hx function among individual fishes at the time of sampling is unclear. Regardless of mRNA abundance on northern blot, full-length Hx cDNA was readily amplified from liver RNA for all icefish species (Fig. 2), indicating that this heme scavenging gene is actively transcribed regardless of presence or absence of Hb. Although Hb would be the largest reservoir of heme, in its absence in icefishes, demands for Hx function could arise from a variety of other cellular client hemoproteins. These include cytochromes of the respiratory chain in mitochondria, which occur at greatly increased density in icefish cells than in the red-blooded species to facilitate cellular O₂ transport, heart Mb in icefish species that synthesize it (O'Brien and Mueller, 2010), microsomal and mitochondrial cytochrome P450s, as well as neuroglobin (Cheng, et al., 2009) and cytoglobin (Cuypers, et al., 2017). This level of client demand appears sufficient to maintain steady levels of Hx transcription. The Hx cDNA sequences of all the icefish species do not contain any deleterious mutations and are thus likely translated into functional Hx proteins (Fig. 2). Our sequence and structural analyses support the functionality of the encoded Hx protein. The strong conservation of the distinctive dual Hx domains with the sequence repeats and functional residues (Fig. 2, Fig. S4) and of the global and domain structures with known Hx (Fig. 7) across red-blooded notothenioids and the icefishes serve to support that proteins encoded by their Hx cDNA would be functional. The difference in expression of the Hp and Hx genes therefore lies in the fact that only the Hp function experiences truly relaxed selection, thus drastic trait reduction was only observed for Hp.

In contrast to Hp, and in keeping with Hx genic integrity and transcription, estimates of the strength of selective pressure showed an intensification in both the icefishes and red-blooded Antarctic notothenioids relative to the reference teleost species (Table 2). Heme occurs in four (a, b, c and d1) common types in vertebrate hemoproteins, with heme b being the most abundant as the prosthetic group of Hb and Mb (Chapman et al., 1997). The detected intensification in selective pressure among the Antarctic notothenioids may again reflect the absence or reduced incidence of free heme from Hb in the white-blooded and red-blooded species, respectively, and the resulting optimization of the Hx protein in the face of what is likely a marked change in the prevalence of heme types.

Alternatively, this could represent a shift in the gene primarily responsible for heme binding in the Antarctic notothenioids. Most teleost species have two WAP65 (Warm Acclimation Protein 65) genes that are similar to the mammalian Hx. Although differences in their physiological function remain unclear, both retain heme binding motifs (Machado et al., 2014; Diaz-Rosales et al., 2014). Only one of these could be found in the available transcriptomes and genomes of the Antarctic notothenioids (WAP65-2), suggesting that the other has been lost. The detected signature of intensified selective pressure could thus reflect this gene now being co-opted into the additional roles formerly carried out by WAP65-1.

Underlying Hp expression loss is a more complex and interesting pattern of molecular lesions that are independent of phylogenetic relatedness, as summarized in Fig. 8 from cDNA (Fig. 1) and genomic sequences (Fig. 3). Three distinct nonsense mutations (at amino acid positions 24, 133 and 174) occurred in five species – P. macropterus, P. georgianus, D. hunteri, C. aceratus and C. antarcticus – occupying separate branches in the icefish phylogeny (Fig. 8). The D. hunteri Hp gene additionally sustained a 7-nt deletion (at positions 52–54), causing a reading frame shift and truncation of the CDS shortly downstream. Any of these mutations will lead to a truncated protein if translated, compromising or obliterating the Hp function. Thus, non-functionalization of the Hp trait upon Hb loss has involved multiple independent molecular mechanisms akin to the independent loss of Mb in various icefish lineages.

Interestingly, despite CDS mutation(s), the Hp gene of P. macropterus, P. georgianus, D. hunteri and C. aceratus continued to be transcribed (Fig. 1), perhaps owing to the stochastic nature of the degeneration process. The combination of intact and pseudogenized Hp genes among the icefishes suggest we are observing differing progression of non-functionalization among these species. The Hp gene, and perhaps others like it, were not lost with the initial disappearance of Hb, but rather their sequences are drifting towards functional extinction in a sporadic manner, and independently among the icefish lineages. Hp is synthesized as a preproprotein, with an N-terminal signal peptide that targets it to the secretory pathway and extracellular export, and proteolytic cleavage of the proprotein that generates the biologically active mature protein (Wicher and Fries, 2006). An interesting observation of note is that the nonsense mutation at amino acid position 24 in C. aceratus Hp is located in the pro sequence of the preproprotein. The predicted proprotein cleavage sequence is ‘RSRR’ at sites 38–41, with the Met on 42 being the start of the mature protein (Fig. 1) (Wicher and Fries, 2006). Thus, the entire mature protein CDS remains free of deleterious mutations, and in principle, the Hp gene of C. aceratus could evolve to encode a cytoplasmic protein should there be selection pressure for its continuance.

Apart from changes in the protein coding region, functional diminution of Hp could also have resulted from changes in the regulatory elements that control gene expression. Obtaining the Hp genomic sequences from large DNA insert BAC clones for the basal red-blooded E. maclovinus as an ancestral proxy – the red-blooded Antarctic D. mawsoni and the icefish C. aceratus – allowed us to investigate the evolution of the control region, and the impact of relaxed selective pressure on cis regulatory elements controlling Hp expression. The core promotor, the TATA box, remains intact in the icefish and identical to that of the two red-blooded species (Fig. 5). Several transcriptional factors (TFs) are known to have essential roles in affecting liver-specific gene expression and in controlling the activity of Hp as part of the acute phase response (APR) in inflammation. These include HNF1, HNF3, HNF4, HNF6, C/EBP and DBP (Pelletier et al., 1998; Schrem et al., 2002; Odom et al., 2004). Among these, the induction of Hp during the APR is itself regulated through C/EBP beta and delta in hepatocytes (Baumann and Gauldie, 1994).

The prolonged loss of selective pressure for Hp in icefishes would expectedly provide the opportunity for sequence degradation in the cis regulatory elements that these TFs interact with, leading to the drastic reduction of Hp expression that we observed. Surprisingly, HNF1A, HNF4A, CEBP-B, CEBP-D and DBP remained conserved in C. aceratus (Fig. 5). Persistence of these elements, however, does not necessarily mean that the transcription of the Hp gene is actively engaged by these TFs, as transcription does have multiple levels of control. It only shows that the reduction in selective pressure thus far has not impacted these putative regulatory elements. It remains unclear whether the icefish Hp gene would be transcriptionally responsive during activation of the APR. Regardless, only Hp genes without deleterious mutations in the protein CDS, such as that of C. wilsoni, which was also able to synthesize the protein, would be functionally relevant.

With the partial Hp genomic sequences spanning the 3′ flanking region, we identified an atypical polyadenylation signal sequence specific to the icefish. The A-to-G mutation in the fifth position of the polyadenylation signal sequence in the icefishes (Fig. 4) has been found to result in the most drastic reduction in the accuracy and efficiency of the cleavage of the precursor mRNA and polyA tailing in a systematic examination of the positional effect of point mutations of AATAAA on these processes (Sheets et al., 1990). The altered polyadenylation reduces the pool of polyA+ mRNA, and ultimately reduces the amount of final protein product (Lutz and Moreira, 2011). The 3′ portion of the Hp genomic sequences we obtained for all 16 icefishes in the family showed that five species now share this mutated polyadenylation signal sequence (Fig. 4C). Their appearance in phylogenetically distinct lineages throughout the clade further corroborate relaxation on the necessity to regularly synthesize transcripts, resulting in stochastic loss.

Many genes play partnering roles in supporting Hb and Mb proteostasis. Studies of icefishes in this regard thus far have logically focused on genes in the erythropoiesis pathways, in part to identify mechanisms behind the erythrocyte non-production in these fish. An early study examined one candidate gene, Bty (blood thirsty) (Detrich and Yergeau, 2004). In the zebrafish test system, Bty protein expression was found to be required for erythrocyte production and Hb synthesis (Yergeau et al., 2005). This was expressed at much lower levels in the icefish C. aceratus than in the red-blooded N. coriiceps; however, the integrity of Bty in icefishes and the basis of its paltry expression remained unknown.

More recent work on the evolution of the erythropoiesis network investigated evolution of transcriptional activity in hematopoietic tissues in the icefish C. hamatus (Xu et al., 2015). This study identified substantial remodeling of the hematopoietic programs in icefish with large-scale reduced expression of hematopoiesis-related genes compared with a red-blooded notothenioid, including many transcriptional factors essential for erythropoiesis. However, in investigating two icefish species, Desvignes et al. (2016) found broad conservation of miRNAs that regulate erythropoiesis in white-blooded Antarctic icefish, despite their lacking RBCs, and suggested possible additional roles outside of erythropoiesis.

In these studies, the focus on single species or small numbers of species limits understanding the evolution of these systems in the diversification of the icefishes, perhaps missing broader evolutionary trends. In addition, no study thus far has considered the fate of those genes associated with the aftermath of RBC senescence, the recycling of otherwise strongly oxidizing and thus cytotoxic free Hb and heme, which should have been concomitantly impacted by the loss of selective pressures from the disappearance of RBCs and Hb. Therefore, the present study is the first to address this hitherto unexplored evolutionary perspective, and assess the functional outcome and the underlying molecular mechanism of two primary supporting genes across the entire icefish family.

The icefishes represent a rich, ongoing ‘evolutionary experiment’. The loss of Hb, the major and abundant oxygen transport protein, simultaneously relaxed selective pressure on the maintenance of suites of supporting genes, and generated novel selective pressures that compel emergence or reprogramming of other suites of genes to sustain the unique Hb-less physiology and life history. The fate of the two relevant genes Hp and Hx contributes to the understanding of gene fate under relaxed selection in the Hb-less state. Although Hx has lost its primary client protein, remaining client proteins appear to exert sufficient selective pressure to ensure Hx functional persistence. In contrast, it is unknown whether teleost Hp has non-Hb client proteins. The Hp gene in several icefishes has sustained deleterious mutations, and the protein is undetectable for all except one tested species, which strongly suggest that should it play any accessory role it is either minimally essential or readily replaceable. Thus, Hb is very likely the dominant, if not the sole client, and its loss would have all but fully relaxed selective pressure on the maintenance of Hp. This would allow the gene to drift down various mutational paths, ultimately leading to the widespread expression loss and non-functionalization observed across the icefish family, exemplifying a clear case of co-evolutionary loss caught in action. This work shows that a more complete understanding of the impacts of Hb loss requires investigations into the dual sides of this impact, both the creation of new selective pressures for compensatory adaptations and the relaxation of selective pressures on former supporting genes. Further, the complex pattern of loss seen among the icefishes suggests that the impacts of relaxed selective pressure remain ongoing, underscoring the importance of family-wide assessments of evolutionary impacts to achieve comprehensive understanding.

We would like to thank the staff and contractors at Palmer Stations and aboard the RV LM Gould for their assistance in carrying out this project. We would further like to acknowledge the work of two undergraduate students, Elizabeth Kalmanek and Nicholas Bart, which contributed to preliminary amplifications of Hp and Hx sequences. Finally, the authors would like to thank Dr Daniel Macqueen for thoughtful comments and critiques that helped refine this work.