Phenotypic flexibility in respiratory traits is associated with improved aerial respiration in an amphibious fish out of water

The transition from an aquatic to terrestrial environment imposes many respiratory challenges for amphibious fishes (Brown et al., 1992; Sayer, 2005). As a result, amphibious fishes have evolved specific adaptations for life out of water (Graham, 1997). Many amphibious fishes switch their primary site of O₂ uptake from the gills to air-breathing organs (e.g. gas bladder, buccal–pharyngeal cavity, skin) to maintain O₂ demands (Graham, 1997). Epidermal capillaries have also been observed close to the skin surface (1–119 µm) in amphibious fishes (Mittal and Munshi, 1971; Grizzle and Thiyagarajah, 1987; Park et al., 2006), whereas in most fishes, capillaries are located deeper within the dermis (Feder and Burggren, 1985). The description of morphological adaptations for air breathing has a long history (e.g. Das, 1934; Hughes and Munshi, 1968), but less attention has been focused on reversible phenotypic flexibility in fishes out of water (Wright and Turko, 2016).

There is some evidence that amphibious fishes enhance aerial respiration out of water (emersion) by altering the efficiency of O₂ uptake. Aerial respiration may require some fish to undergo structural modifications (e.g. reduction in the diffusion distance or an increase in the number of cutaneous capillaries, i.e. angiogenesis) to maximize surface area for exchange (Marusic et al., 1981; Cooper et al., 2012; Glover et al., 2013; Turko et al., 2014). Angiogenesis is the development of new capillaries derived from pre-existing blood vessels (Djonov et al., 2000). This process can be controlled via different mechanisms (i.e. capillary intussusception and sprouting) and through multiple genes [i.e. vascular endothelial growth factor (vegf), angiopoietin-1 (angpt1), angiopoietin-2 (angpt2), ephrins (efn); Prior et al., 2004; Fagiani and Christofori, 2013]. The gene coding for platelet endothelial cell adhesion molecule (pecam-1) expresses a protein (CD31) that helps form junctions between endothelial cells (Albelda et al., 1991). Angiogenesis during emersion would presumably increase blood flow near the respiratory epithelium, maximizing gaseous exchange (Glover et al., 2013).

Plasticity in O₂ transport may also play a role in enhancing respiration in amphibious fishes out of water. For example, a faster rate of blood delivery (increased heart rate) would increase O₂ transport, as reported in mudskippers, Periopthalmodon australis (Kok et al., 1998; Garey, 1962). Reversible phenotypic plasticity of haemoglobin (Hb) properties would also ameliorate the impact of CO₂ accumulation and blood acidosis in amphibious fishes out of water. Increased Hb–O₂ affinity during air exposure may be beneficial in offsetting the Bohr shift owing to CO₂ retention in emersed amphibious fishes (Graham, 1997; Morris and Bridges, 1994). By altering O₂ carrying capacity (Hb concentration and/or erythrocyte density; Delaney et al., 1976; Johansen et al., 1976; Marusic et al., 1981; Urbina and Glover, 2012; Turko et al., 2014), amphibious fishes may compensate for reduced O₂ carrying capacity (Root effect; Root, 1931). Some fishes do both. For example, Kryptolebias marmoratus, Protopterus amphibious and Protopterus aethiopicus increase their Hb–O₂ affinity (lower P₅₀) and increase Hb concentration during emersion (Delaney et al., 1976; Johansen et al., 1976; Turko et al., 2014). Thus, rapid responses to enhance O₂ transport would offset the negative effects of elevated blood CO₂ in air-exposed fishes (Graham, 1997).

We tested the hypothesis that amphibious fishes that leave water, have enhanced respiratory performance on land as a result of rapid functional phenotypic flexibility of respiratory traits. This hypothesis predicted that fish with increased cutaneous angiogenesis in response to air exposure would have a higher terrestrial respiratory performance [lower critical oxygen tension (P_crit), higher regulation index (RI)]. As well, fish that have increased blood carrying capacity [increased hematocrit (Hct), increased number of red blood cells (n_RBC)] in response to air exposure should have a higher terrestrial respiratory performance (lower P_crit, higher RI). P_crit is defined as the point at which the O₂ consumption rate of an organism becomes dependent on environmental O₂ levels (Ultsch et al., 1978). In aquatic environments, P_crit has been found to be highly correlated to respiratory traits along the O₂ cascade (e.g. gill surface area and Hb–O₂ affinity; Mandic et al., 2009). In contrast, RI, an alternate performance measure, is the overall regulatory ability of the fish over the full range of atmospheric O₂ levels. This parameter provides insight as to whether an organism is more a conformer or a regulator (Mueller and Seymour, 2011).

Kryptolebias marmoratus is an ideal species for studying the terrestrial respiratory performance of amphibious fish because they can tolerate weeks out of water (Taylor, 2012; Taylor et al., 2008; Wright, 2012), and prolonged air exposure results in angiogenesis of alternate respiratory surfaces (Cooper et al., 2012; Turko et al., 2014) and increased blood Hb concentration (Turko et al., 2014). In addition, K. marmoratus are one of only two known self-fertilizing vertebrates, creating isogenic offspring (Harrington, 1961), which allowed us to control genetic variation while manipulating the environment. Therefore, we compared respiratory traits (Hct, n_RBC and angiogenesis) and performance (aerial P_crit and RI) across multiple isogenic lineages of K. marmoratus.

Kryptolebias marmoratus Poey 1880 hermaphroditic strains were obtained from the breeding colony housed at the Hagen Aqualab at the University of Guelph, Guelph, Ontario, Canada. The isogenic strains of FW2 (freshwater) (Platek et al., 2017), 50.91 (Belize), HON11 (Honduras) and SLC (Florida) were used (Tatarenkov et al., 2010). Fish were held individually in 120 ml semi-transparent plastic containers (FisherBrand Collection Containers, Fisher Scientific) and maintained under constant conditions [12 h:12 h light:dark cycle, 25°C, 15 ppt salinity for Belize, Honduras and Florida strains (Frick and Wright, 2002) and 0.3 ppt salinity for freshwater fish (Platek et al., 2017)]. Brackish water and freshwater were made with reverse osmosis water and marine salt (Instant Ocean, Crystal Sea) to the appropriate salinity and changed weekly. The fish were fed Artemia nauplii three times a week until the beginning of experiments. This project was approved by the University of Guelph Animal Care Committee (AUP 2239).

Fish [Honduras (0.120±0.003 g), Belize (0.120±0.003 g), Florida (0.120±0.003 g) and freshwater (0.140±0.004 g)] were acclimated to water (control) or air (1, 3 or 7 days) at 25°C. Air-acclimated fish were maintained on moist filter paper (15 or 0.3 ppt) in plastic containers, as previously described (Ong et al., 2007). All treatment groups were subjected to a P_crit test and an RI test in air in order to determine their respiratory performance ability. It was necessary to perform all experiments under the same conditions for comparisons. It is likely that some changes occurred very quickly when K. marmoratus left the water (1–2 h), but these potential changes are presumably minor relative to the more profound changes observed at 1, 3 and 7 days. Owing to the small size of the fish, all measurements could not be performed on the same individuals. New groups of fish were acclimated to air or water as described above and used for histological, blood or gene expression analyses (see below and Table S1). At the end of the experiment, fish were euthanized with tricaine methanesulfonate (MS-222; 1.5 mg ml⁻¹) and cut into transverse sections anterior to the dorsal fin. Sections were covered with embedding medium (Shandon cryomatrix, ThermoFisher Scientific), frozen using liquid-nitrogen-chilled 2-methylbutane and stored at −80°C until sectioning (Brunt et al., 2016). An additional experiment on separate fish was conducted to measure Hct and n_RBC. Fish from each strain were air exposed for 0, 1 or 7 days. Blood was collected by caudal severance using heparinized microhematocrit tubes (Kimble Chase) (Turko et al., 2014). For gene expression analyses, fish (Honduras and Florida strains only) were sampled at 0 h (pre-emersion) and post-emersion at 1 h, 6 h, 1 day, 3 days and 7 days. Skin samples were immediately transferred into RNAlater (ThermoFisher Scientific) and stored at −20°C until analyses. It is important to note that we chose to only analyze the Honduras and Florida strains for gene expression analysis based on their differences in survival rate out of water. Honduras fish had a significantly higher survival rate after 7 days compared with freshwater fish (Y.W.D., T.S.B., J. Schmitz, S. Kelly, P.A.W. and A.W., in preparation). Additionally, the cost of RNA-seq limited the number of strains we could analyze.

P_crit was measured in custom-made glass micro-respirometry chambers (∼1 ml) in which an optode was used to measure O₂ saturation (Loligo Systems WITROX 4). Chambers were kept in an incubator (Innova 4230, New Brunswick Scientific) to maintain a constant temperature of 25°C. Before each experiment, wet filter paper was inserted into each respirometry chamber to maintain a humid environment during air exposure. All experiments were conducted between 12:00 and 18:00 h to account for diurnal fluctuations in metabolic rate (Rodela and Wright, 2006). Preliminary experiments were conducted to determine the appropriate ratio of mass of fish to volume of chamber to achieve a significant change in atmospheric O₂ in a reasonable period of time. The volume of the chamber was adjusted by adding an inert material (wax). P_crit in air was measured using a modification of a step-wise hypoxia protocol previously described for an aquatic system (Borowiec et al., 2015; Crans et al., 2015), with a few exceptions to account for the differences in O₂ content between water and air. Fish were inserted into respirometry chambers and were acclimated to the chamber for 20 min at 100% air saturation. Preliminary experiments showed that 20 min was a sufficient acclimation period in air. This was determined by measuring the rate of O₂ consumption (Ṁ_O2) over a 2-h period, and we found that there was no statistically significant decrease or change across all Ṁ_O2 time points. Additionally, in preliminary experiments conducted on measuring maximum metabolic rate out of water, we found that metabolic rate decreased back to resting metabolic rate in only a few minutes (T.S.B., unpublished data), which supports the idea that K. marmoratus recover quickly in response to handling stress.

O₂ consumption was initially measured at a partial pressure of O₂ (P_O2) of 21.2 kPa and then at a P_O2 of 14.8 to 10.6 kPa in steps of 2.1 kPa and from 10.6 to 1.1 kPa in steps of 1.1 kPa. At each step, the difference in P_O2 was recorded over 10 min in a sealed chamber. After each measurement, the chamber was flushed with the new P_O2 air and fish were left for 5 min before the next measurement began. Control of the P_O2 was achieved using a gas mixing system with air and N₂ (Wosthoff, Calibrated Instruments Inc.). Optodes were calibrated weekly using air (100% P_O2) and 2 mol l⁻¹ of sodium sulfite (0% dissolved O₂) as described previously (Sutton et al., 2018). Routine metabolic rate (RMR) was calculated as μmol O₂ g⁻¹ h⁻¹ by measuring the slope of the O₂ consumption curve over time at 21.2 kPa. Background respiration was measured before and after each experiment; however, it was found to be negligible.

Aerial P_crit for each fish was calculated using nonlinear regressions, as described by Marshall et al. (2013), which better accommodates data sets in which RMR more gradually declines with environmental O₂ levels, rather than a sharp transition at a specific O₂ level. To calculate RI, we first determined the curve that best fit the data and then we fitted a straight line at the start and end of Ṁ_O2. From there, we calculated RI as the area between the curve and straight line as described by Mueller and Seymour (2011). An RI of 1 represented complete regulation and a value of 0 represented total conformity to environmental O₂ levels (Mueller and Seymour, 2011).

Immunofluorescence was used to stain for the endothelial intercellular junction protein cluster of differentiation 31 (CD31), as previously described (Cooper et al., 2012), with a few modifications. The CD31 antibody has previously been used in the literature to quantify changes in endothelial cell proliferation in fishes (Cao et al., 2008; Cooper et al., 2012) and angiogenesis in mammals (DeLisser et al., 1997). Frozen transverse sections were cut (8 µm thick) using a cryostat at −22°C (Leica CM3050 S) and slides were stored at −80°C until staining. Slides were defrosted for 2 h prior to staining and then rinsed in phosphate-buffered saline (PBS) with Triton-X (0.1% v/v) for two 5-min washes to permeabilize the tissues. Samples were blocked for 1 h at room temperature in blocking solution [PBS, 5% normal goat serum, 0.1% (v/v) Tween-20, 0.05% (v/v) sodium azide]. All samples were incubated in a humidified chamber overnight at 4°C in primary antibody [1:100 rat anti-mouse PECAM/CD31:PBS (cat. no. 553370, BD Pharmingen), 0.1% (v/v) Tween-20, 0.05% (v/v) sodium azide]. Samples were then rinsed in PBS with Tween-20 (0.1% v/v) three times for 5 min each. Samples were incubated in a humidified chamber for 2 h at room temperature with Alexa-Fluor-488-labeled secondary antibody (1:400 goat anti-rat IgG:PBS; Invitrogen). Samples were washed five times for 5 min each with PBS and mounted with Fluoromount with DAPI (Sigma-Aldrich). A negative control, in which no primary antibody was applied, was used to ensure the specificity of the secondary antibody. Images were taken (20×) on the same day using a Nikon epifluorescent microscope (Nikon Eclipse 90i microscope) using the same camera settings for all images. Using ImageJ, a line was traced around the epidermis (dorsal or ventral) and the integrated fluorescence density was calculated. To account for potential differences in tissue thickness across samples, values were normalized by dividing the integrated density of the epidermis by the integrated density of the skeletal muscle.

Immediately following skin dissections, tissues were preserved in RNAlater and archived at −20°C. Skin tissues (N=5 per time point) were individually homogenized (Omni BeadRuptor) in RLT buffer with β-mercaptoethanol and RNA was purified from homogenates using Qiagen RNeasy Purification kits. RNA-seq libraries were prepared using NEBNext RNA library preparation kits for Illumina. Libraries from five replicate individuals per strain, per treatment were prepared. Each sample was tagged with a unique barcode. All samples were multiplexed into a single pool (including samples not analyzed as part of this project), and this pool was sequenced across four lanes of Illumina HiSeq 4000 (PE-150). Sequencing yielded 1,236,473,082 raw reads across the 60 experimental samples. Sequencing failed for two samples, including one sample from the Honduras strain from the 72-h emersion sampling treatment (N=4), and one sample from the Honduras time-0 immersion control treatment (N=4). Short and low-quality reads were removed with Trimmomatic 0.36 (Bolger et al., 2014). Reads were mapped to the reference genome (RefSeq assembly accession: GCF_001649575.1) using STAR (Dobin et al., 2013). The average number of mapped reads per sample was 20,414,142. Read counts were generated using HTSeq (Anders et al., 2015). We removed genes from subsequent analyses when read counts were too low (criteria: read counts >10 in at least 5 samples). Read counts were log₂ transformed and normalized for gene length and total library size in edgeR (Robinson et al., 2010).

Hct was measured after centrifugation (International Clinical Centrifuge, Model CL, International Equipment) at 5200 g for 2 min. Because blood volumes were minute (<1 µl per fish), images were taken of the microhematocrit tubes using a dissecting microscope (Wild of Canada Limited) and the proportion of packed red blood cells was determined using ImageJ (Bianchini and Wright, 2013). To measure the number of red blood cells (n_RBC), whole blood was diluted in a 1:400 dilution (whole blood:Cortland's isotonic saline) (Wolf, 1963) and then further diluted 1:1 with 0.4% Trypan Blue solution to stain for non-viable red blood cells (Turko et al., 2014). Red blood cells were counted using a standard hemocytometer (American Optical) using a Nikon Eclipse 90i epifluorescent microscope. Red blood cells were manually counted from a single row in the center square of the hemocytometer. Rows were randomly selected by assigning each row a number and using a random number generator to determine which row to count. Unfortunately, we were unable to measure Hb–O₂ affinity as in previous studies (Bianchini and Wright, 2013; Turko et al., 2014) owing to equipment failure.

The aerial O₂ consumption rate of all four strains decreased with decreasing atmospheric P_O2 (Fig. 1). Time and strain had significant effects on critical O₂ tension measured as aerial P_crit (time: F_3,156=34.13, P<0.001; strain: F_3,156=3.14, P=0.02; Fig. 2A), but their interaction was not significant (time×strain: F_9,156=1.53, P=0.14). Aerial P_crit was significantly lower at 1, 3 and 7 days of air exposure relative to control (0 days), and was also significantly lower at 7 days relative to 3 days of air exposure. The Belize strain had a significantly higher aerial P_crit than the Honduras strain (P<0.01), whereas the Florida and freshwater strains were intermediate (P>0.05; Fig. 2A).

Air exposure had a significant effect on RI. Fish acclimated to air for 1 day had a significantly higher RI than fish acclimated to air for 0 and 7 days, whereas fish acclimated for 3 days had an intermediate RI (time: F_3,156=3.90, P=0.01; Fig. 2B). There was no effect of strain on RI (strain: F_3,156=0.22, P=0.88).

RMR was also altered in fish out of water. RMR was significantly lower after 7 days in air relative to 0, 1 and 3 days of air exposure at 21.2 kPa, and this was not influenced by strain or their interaction (time: F_3,156=7.62, P<0.001; strain: F_3,156=0.47, P=0.70; interaction: F_3,156=1.54, P=0.14; Fig. 2C).

Angiogenesis was enhanced by air exposure. CD31 expression was visible in both the dorsal and ventral region of the epidermis. Moreover, the expression of CD31 appeared more prominent at 3 and 7 days of air exposure in both regions (Fig. 3A–H). In the dorsal region of the epidermis, fish acclimated to air for 1 day had a significantly higher CD31 fluorescence intensity than fish acclimated to air for 0 days. Moreover, fish acclimated to air for 3 and 7 days had a significantly higher fluorescence intensity relative to fish acclimated to air for 0 and 1 day (time: F_3,105=22.74, P<0.001; strain: F_3,105=0.15, P=0.93; interaction: F_9,105=0.73, P=0.68; Fig. 3I). In the ventral region of the epidermis, fish acclimated to air for 3 and 7 days had a significantly higher fluorescence intensity than fish acclimated to air for 0 and 1 day (time: F_3,107=7.29, P<0.001; strain: F_3,107=0.70, P=0.55; interaction: F_9,158=0.33, P=0.96; Fig. 3J).

Angiogenesis was linked to respiratory performance across strains. The change in dorsal angiogenesis was positively related to the change in RI at 1 and 3 days of air exposure (both P=0.01); however, no significant relationship was detected at 7 days of air exposure (Table 1).

We found a strong upregulation of three angiogenesis genes in the skin after 1 h in air (vegfa 1.9-fold; angpt2 3.7-fold; efna 7-fold; P<0.05; Fig. 4A–C) compared with the control. In contrast, pecam-1 was significantly upregulated (1.6-fold; P<0.05; Fig. 4D) by 6 h following emersion. However, by 7 days of air exposure there was no significant difference in expression across all angiogenesis genes compared with the control.

Hct was altered by air exposure. Both strain and air exposure time had direct and interacting effects on Hct (time: F_2,82=3.28, P=0.04; strain: F_3,82=3.73, P=0.01; interaction: F_6,82=4.27, P<0.001; Fig. 5A). At day 0, the Belize strain had a significantly higher Hct relative to both the Florida and Honduras strains, and the Hct of the freshwater strain was also higher than that of the Florida strain. Only the Florida strain showed a significant increase in Hct after 1 day of air exposure compared with 0 days, but both the Honduras and Florida strains had higher Hct values after 7 days. Neither air exposure nor strain influenced the number of red blood cells (all P>0.05, Fig. 5B).

Initial changes in respiratory performance were related to O₂ carrying capacity. The change in aerial P_crit was positively and significantly related to the change in Hct at 1 day of air exposure (P=0.01; Table 1). No other significant relationships were detected between O₂ carrying capacity (Hct, n_RBC) and respiratory performance variables (aerial P_crit, RI; Table 1).

In this study, we experimentally demonstrated that aerial acclimation improves respiratory performance in an amphibious fish. We compared respiratory performance out of water (aerial P_crit, RI) with traits linked to O₂ uptake (cutaneous angiogenesis) and O₂ transport (Hct, n_RBC) in four isogenic lineages of K. marmoratus. In general, we found that respiratory performance and traits varied across both emersion time and strain, but time was the stronger factor. Indeed, we found that time had the largest effect size in five of the seven traits measured, except for the two blood parameters (Table 2). Kryptolebias marmoratus showed a consistently improved respiratory performance (lower aerial P_crit and higher RI) after only 1 day of aerial acclimation. The initial rapid improvement in aerial P_crit was most strongly linked to O₂ carrying capacity (Hct), and the initial improvement in RI was significantly associated with dorsal angiogenesis. These results suggest that aerial P_crit and RI are regulated by different factors along the O₂ cascade. Overall, K. marmoratus showed modifications along both the O₂ uptake and transport system in response to air over time; however, only initial plastic changes in respiratory traits were related to improved respiratory performance out of water.

There appeared to be two phases of improved respiratory performance occurring over time in air. An initial rapid response (1 day; lower aerial P_crit and higher RI) was followed by a delayed response (7 days; lower aerial P_crit), suggesting that two different physiological responses may be involved. Both aerial P_crit (lower) and RI (higher) were significantly improved by 1 day in air. Interestingly, aerial P_crit stabilized between 1 and 3 days of emersion but was even lower by day 7. In contrast, RI returned to control values by day 7. Moreover, the Honduras strain had an overall different aerial P_crit than the Belize strain across all time points; however, the largest variation in aerial P_crit appears to be at day 1. The Honduras strain appeared to have a more rapid respiratory response relative to the Belize strain. Interestingly, a previous study in our laboratory found that the K. marmoratus Honduras strain survived significantly longer in air compared to both the Florida and Belize strains (A. Turko, J. Doherty, P. A. Wright, unpublished data). Therefore, a higher emersion tolerance may be related, in part, to a greater initial respiratory ability during emersion.

The data indicate that the initial rapid plastic change in aerial P_crit was primarily driven by O₂ carrying capacity. Our hypothesis predicted that fish with increased O₂ carrying capacity would have a higher respiratory performance (lower aerial P_crit and higher RI). In support of this, we found that the change in Hct was strongly correlated (R²=0.97) to the change in aerial P_crit at 1 day of air exposure, where strains with the largest increase in Hct also had the largest decrease in aerial P_crit. The relationship between Hct and P_crit has been shown in other species; for example, hypoxia-tolerant aquatic fishes with a higher Hct tend to have a lower P_crit (Chapman et al., 2002). However, Mandic et al. (2009) found no significant relationship between Hct and P_crit in various species of marine sculpins. Thus, the importance of O₂ carrying capacity in fish respiration may be species, time or environment dependent.

The observed initial relationship between Hct and aerial P_crit in K. marmoratus may be influenced by inherent differences in Hct across strains or the ability to modify Hct in response to air exposure. Baseline O₂ carrying capacity was significantly different across strains: both the Belize and freshwater strains had higher Hct relative to the Florida strain. Elevated Hct is thought to be important for both O₂ uptake (gills or skin) and O₂ delivery to the tissues (Wells et al., 2003), as well as mitigating the effects of elevated tissue CO₂ on O₂ carrying capacity (Graham, 1997). Thus, having an inherently high O₂ carrying capacity may be beneficial in sustaining O₂ demands during the initial transition onto land if the onset of respiratory plastic changes is delayed. Additionally, strains that had a higher baseline O₂ carrying capacity exhibited no change in Hct in response to aerial acclimation. Increased Hct can be also a disadvantage as there is an exponential increase in blood viscosity with small changes in Hct, which can hinder blood flow through the epidermal capillaries and increase work output by the heart (Baldwin and Wells, 1990; Wells and Weber, 1991). Therefore, the costs associated with a higher Hct may exceed the advantages of having a higher O₂ carrying capacity in strains with elevated Hct (Urbina and Glover, 2012).

The mechanisms involved in the rapid increase in Hct observed during emersion are unknown. Acute changes in Hct can arise as a result of a shift in fluid volume or through the release of red blood cells via the spleen (Jensen et al., 1993; Gallaugher and Farrell, 1998). An acute increase in Hct after 1 day of air exposure was only observed in the Florida strain, but there was no significant change in the n_RBC. The increase in Hct could reflect a change in cell volume (Weber and Jensen, 1988). We estimated mean cell volume (MCV) from dividing mean Hct by mean n_RBC (Turko et al., 2014). Indeed, an ∼11% increase MCV was found in the Florida strain, whereas in the other strains, MCV tended to decrease (Honduras, Belize) or the change in cell volume was negligible (∼1% for freshwater; data not shown). Whether elevated catecholamines were involved in the Hct changes in the Florida strain is unknown, but β-adrenergic-stimulated erythrocyte swelling can result in response to air exposure (Nikinmaa, 1982; Perry et al., 1989). Further investigation will be necessary to tease apart the mechanism involved.

An increase in the number of blood vessels and/or epidermal blood perfusion may improve the ability to regulate O₂ consumption during emersion. We predicted that strains with higher cutaneous vascularization would exhibit a higher respiratory performance during emersion. Indeed, there was a strong relationship (R²=0.97) between the change in dorsal angiogenesis and RI at 1 and 3 days of air exposure, where strains that showed the largest increase in CD31 expression also showed the largest increase in RI. In a separate study on K. marmoratus, RI was associated with gill surface area, implying that the ability to regulate O₂ consumption (RI) may be linked to modifications across the O₂-uptake system (Turko et al., 2012). This rapid adjustment in the dorsal region of the epidermis is also consistent with recent behavioral data showing that K. marmoratus spent significantly more time exposing their dorsal surface to air relative to their ventral or lateral sides (Heffell et al., 2017). Mudskippers (Boleophthalmus and Scartelaos) also have a high degree of vascularization in the head and dorsal area, regions most often exposed to air (Zhang et al., 2000).

The rapid changes we observed in angiogenesis were consistent with the dramatic increase in the expression of genes involved with blood vessel development. Angiogenesis requires complex multi-step signaling that is mediated by molecules from three protein families – vascular endothelial growth factors (VEGFs), angiopoietins (ANGPTs) and ephrins (EFNs) – that act through receptor tyrosine kinases in endothelial cells (Gale and Yancopoulos, 1999). Transcripts for key members of each of these ligand families were significantly upregulated within 1 h of air exposure in K. marmoratus skin (Fig. 4). Parallel upregulation of both angpt2 and vegfa is consistent with promotion of angiogenesis (Maisonpierre et al., 1997; Holash et al., 1999). VEGF-induced angiogenesis requires ephrin-A2 (EPHA2) receptor activation, and VEGF induces expression of ephrinA1 (the ligand for EPHA2) in endothelial cells (Cheng et al., 2002). We found that 7-fold upregulation of epha2 occurs within 1 h of emersion. The gene that codes for CD31, pecam-1, was increased in expression by 1.6-fold within 6 h of emersion. This delay relative to the other angiogenesis-associated genes may reflect the time required to grow new endothelial cells before cell-to-cell junctions are formed. It is important to note that CD31 (pecam-1) is also present in lymphocytes, platelets, leukocytes (neutrophils) and monocytes; thus other physiological roles of CD31 include involvement in the inflammatory response and vasculogenesis during embryonic development (DeLisser et al., 1994; Pinter et al., 1997). However, it is unlikely that our data are signaling these other physiological processes. Together, these data indicate that key initiators of angiogenesis signaling pathways are coordinately upregulated almost immediately following exposure to air in K. marmoratus skin. Moreover, the gene expression and immunofluorescence data are strong evidence that K. marmoratus exhibit cutaneous angiogenesis during emersion, potentially as a mechanism to increase cutaneous gas exchange as well as the transfer of other molecules.

Respiratory traits were not significantly linked to variation in either aerial P_crit or RI after 7 days in air. In fact, variation in cutaneous angiogenesis and P_crit were low between strains, although both were significantly enhanced relative to earlier time points across all strains. Cutaneous angiogenesis may be more important for other functional mechanisms (i.e. ion, water and nitrogen regulation) rather than respiration during more prolonged emersion. It is also possible that cutaneous respiration was augmented by angiogenesis in the buccal/opercular regions, as K. marmoratus are known to occasionally gulp air out of water (Turko et al., 2014). Finally, other physiological factors within the O₂ transport system may also play a role in the improved respiratory performance after 7 days in air (e.g. Hb–O₂ affinity; Johansen et al., 1976; Turko et al., 2014).

The O₂-transport system (Hct) may be less important in long-term improvement of respiratory performance during emersion. An increase in O₂ carrying capacity was found in both the Honduras and Florida strains, where both strains had a significantly higher Hct at 7 days of air exposure relative to the control values, but no change in n_RBC was observed. Thus, increased Hct could be the result of a shift in plasma volume (Gallaugher and Farrell, 1998). In contrast to the findings at 1 day of air exposure, there was no significant relationship between the change in aerial P_crit and Hct at 7 days of air exposure, nor was there a significant relationship between RI and the change in Hct or n_RBC.

RMR in response to air exposure is highly variable across air-breathing fish species. Some species increased (Gordon et al., 1970; Sacca and Burggren, 1982; Urbina et al., 2014), decreased (Delaney et al., 1974; Berg and Steen, 1965; Tamura et al., 1976) or did not change (Gordon et al., 1969; Pelster et al., 1988) aerial O₂ uptake compared with aquatic O₂ uptake. At 7 days of air exposure, we found a consistent and significant decrease in RMR, possibly because of a programmed metabolic depression (Storey and Storey, 1990) or the inability to feed during emersion, which would reduce overall energy usage and decrease metabolic demands (O'Connor et al., 2000).

Small sample size in this study limited the statistical power. Owing to the small size of the fish (∼0.12 g) and minute blood volumes, it was not possible to complete all measurements on the same individuals. Although the significant R² values were robust (>0.97), future work should include a larger number of strains and/or other larger amphibious species.

Overall, our findings support the hypothesis that reversible plasticity of the O₂ cascade in amphibious fishes plays a functional role during emersion through the enhancement of respiratory performance. Moreover, we propose that K. marmoratus exhibit two different phases in the enhancement of aerial respiration: an initial rapid response (lower aerial P_crit, higher RI) and a delayed response (lower aerial P_crit and RMR). In turn, these findings may shed some light on the behavior and ecology of K. marmoratus in the wild. Kryptolebias marmoratus display two types of emersion behaviour in the field: (1) short-term emersion, in which they move in and out of water throughout the day to escape either poor water conditions or to capture prey (Taylor, 2012), and (2) long-term emersion that occurs during the dry season, in which they seek moist crevices, including excavated tunnels within decaying mangrove logs for weeks at a time (Taylor et al., 2008). Therefore, the observed initial rapid improvement in respiratory performance would allow K. marmoratus to satisfy metabolic demands during short-term emersions, whereas the delayed reduction in RMR would be an energetic advantage during seasonal long-term fasts out of water.

Aerial P_crit varied across isogenic lineages, suggesting that strains may have different respiratory abilities in air that could affect their emersion tolerance. In turn, this could be the result of differences in their environment or genetics. The largest strain variation we observed was between Honduras and Belize. Thus, we could speculate that the characteristics of their geographic origin may be different (e.g. decreased food availability, longer dry season, fewer predators) and thus might contribute to a more aerial phenotype as displayed in the Honduras strain. However, multiple genetic lineages (as a result of selfing) arise within each population at each geographic site;therefore, it is unlikely that habitat alone contributes to the physiological variation we observed. Self-fertilization is the prevalent form of reproduction in wild populations, as indicated by a high proportion of homozygous individuals, even though outcrossing is also possible (Mackiewicz et al., 2006a, b). Evaluating whether the relative fitness of self-fertilizing lineages could be influenced by respiratory performance in air is a critical next step.

We would like to thank Drs Nick Bernier and Beren Robinson for advice on experimental design and statistical analyses. We would especially like to thank Andy Turko for help in developing the method to measure aerial metabolic rate. We would also like to thank undergraduate student volunteers and work study students, as well as Mike Davies and Matt Cornish for fish husbandry (University of Guelph Hagen Aqualab). We also thank Jennifer Roach for assistance with gene expression data collection (University of California Davis).