Models of primary sex ratios at a major flatback turtle rookery show an anomalous masculinising trend

2.1 Incubation times and embryonic development rate

Of 296 eggs incubated in the laboratory under a range of thermal regimes (see ‘Methods’ section), 156 hatched and almost all eggs that did not hatch showed no signs of development (confirmed by dissection). All 156 hatchlings plus 138 late-stage embryos sampled from field nests were sexed using histological techniques, and the classification of sex was highly repeatable (binomial test; p < 0.0001, probability of the same classification = 0.97).

Incubation times varied between 43 and 60 days (Additional file 1: Figure S1). Eggs incubated predominantly at 32.5°C had a mean incubation duration of 45.4 ± 0.2 days and hatched earlier (unpaired t-test, p < 0.0001) than eggs incubated predominantly at 28°C (57.0 ± 0.4 days). A non-linear development rate function fitted to average incubation times had the parameters b₁ = 2.39, b₂ = 18.03, b₃ = 33.35, b₄ = 6 and b₅ = 0.4 (root mean sum of squares = 3.70; Figure 2).

2.2 Thermosensitive period (TSP)

The controls at 28°C and 32.5°C (treatments 6 and 12) produced 92% male and 100% female hatchlings, respectively (Additional file 2: Figure S1). Female–male–female switches (32.5°C ➔ 28°C ➔ 32.5°C) produced 100% female hatchlings in most treatments (treatments 7, 8, 9 and 11), with the exception of treatment 10 where one of 12 hatchlings was male. Conversely, most male–female–male switches (28°C ➔ 32.5°C ➔ 28°C) produced mixed sexes. Treatments 1 and 3 were switched to female-producing temperatures at 43% of development, thereby inferring a lower boundary of the TSP (Figure 3). Treatment 5, where the switch commenced at 66% of development, produced all male hatchlings (Figure 3), indicating that this switch period fell after the TSP. Hence, the upper boundary of the TSP for a female trigger was 66% of development. As it was not possible to delineate the TSP for a male trigger, the feminising TSP of 43%–66% was hereafter used for estimating primary sex ratios.

2.3 Pivotal temperature (Tpiv) and transitional range of temperatures (TRT)

Using all available sex-ratio data from the field and laboratory, the best fitting TSD function estimated Tpiv at 29.4°C (Figure 4A). When we instead fitted a TSD function to the laboratory data and subsets of field data (see ‘Methods’ section), it produced no differences in the mean residuals (one-way ANOVA; p = 0.997), indicating that fits were effectively the same (Table 2). However, a function fitted to laboratory data and data from a particular subset of field nests (1, 4, 6, 10 and 15) produced the lowest Akaike information criterion (AIC) of 25.9 and the lowest mean residual (Table 2). This function had the parameters P = 29.4, S = −0.02 and K = 0.1, which defines the TRT as between 27.7°C and 31.1°C (Figure 4B). This TSD function was applied for estimating sex ratios following the process outlined in Figure 1.

2.4 Validation of the sand temperature model

Temperature loggers (iButtons) distributed among 11 nests were retrieved from depths between 36 and 80 cm, and some were found clustered together and essentially provided replicate measurements (Additional file 1: Table S1). Thirty-two of 44 loggers deployed were found to contain a complete temperature record (other loggers failed or else recorded temperatures for only part of the incubation period). The diel temperature range decreased with depth, and this was reflected in both the empirical data and the modelled sand temperatures. On average, the difference in the mean measured nest temperatures relative to the mean of the modelled sand temperatures driven by weather data recorded at the rookery ranged from 0.35 to 1.30°C, with an average difference of 0.65°C ± 0.04°C. Sand temperatures modelled using data from the Australian Water Availability Project (AWAP)—a database of historical climate records for Australia [23]—differed from measured nest temperatures by 0.30°C–1.20°C, with an average difference of 0.66°C ± 0.05°C. In general, there was good agreement between actual and modelled sand temperatures (Figure 5D and E) but measured temperatures tended to be slightly lower than NicheMapR estimates near the start of incubation and higher than estimates towards the end of incubation (Figure 5A–C and Additional file 1: Table S1).

2.5 Validation of the physiological model

2.6 Historical sand temperatures and sex ratios

Driving NicheMapR with historical climate data (AWAP: 1990–2013) illustrates the annual variation in sand temperature at a typical nest depth of 50 cm (Figure 6). Annual maximum sand temperatures were stable over this period (linear regression analysis; p = 0.11 and r = 0.34), but minimum sand temperatures decreased (linear regression analysis; p = 0.017 and r = 0.49). Consequently, the annual range of sand temperatures increased between 1990 and 2013 (linear regression analysis; p = 0.004 and r = 0.57). In accordance with these trends, a greater portion of the year favoured the production of males in the most recent years (Figure 6).

Contrasting a cooler year (2011) with a warmer year (2012) in detail (Figure 7) shows that in the former case, mixed sex ratios could be produced over a larger portion of the year. In the warmer year, sand temperatures mostly exceeded the upper boundary of the TRT (31.1°C); hence, females would be produced if the TSP fell within this period. The maximum sand temperature in each year was 34.0°C (December 2011 and November 2012; Figure 7). In both years, sand temperatures were below the lower boundary of the TRT (27.7°C) between May and September, indicating when male hatchlings would be produced.

2.7 Sand temperatures and sex ratios under climate change

Adjusting NicheMapR inputs by increasing the 2007 AWAP air temperatures by the increments projected under 2030 and 2070 emissions scenarios shows that sand temperatures increase by a similar magnitude (Figure 8). Under a 2030 low emissions scenario, sand temperatures increase on average increase by 0.49°C relative to 2007 (Figure 8A). Sand temperatures predicted under the 2030 high emissions scenario and 2070 low emissions scenario are similar, with average increases of 1.22°C and 1.47°C relative to 2007. Under these latter scenarios, conditions that favour male hatchlings are limited to June–September. All other months favour female primary sex ratios (Figure 8B,C). Under a 2070 high emissions scenario, sand temperatures increased on average by 2.78°C, and the only months where temperatures fall below the TRT and thereby favour males are June–August (Figure 8D).

	Sand Depth (cm)
	0–2.5	2.5–10	10–30	30–60
Bulk density (g/cm3)	1.15	1.15	1.24	1.24
Saturated water content (cm3/cm3)	0.4	0.4	0.4	0.4
Mineral thermal conductivity (W/mK)	8.8	8.8	8.8	8.8
Mineral heat capacity (J/KgK)	800	800	800	800
Mineral density (g/cm3)	2.66	2.66	2.66	2.66
% sand moisture	0	6.75	7.57	9.84
Reflectance (%)	41.15	-	-	-

Italics type indicates values measured from sand collected at Cape Domett.

3.1 Validation of sand temperature predictions

The validity of the predicted sand temperatures is strengthened by comparison with actual measurements of nest temperatures at equivalent depths (Additional file 1: Figure S1 and Figure 5). NicheMapR’s microclimate model produced accurate estimates of nest temperatures (without any data-based calibration) across a wide range of depths, except towards the end of incubation when metabolic heating is likely to have increased nest temperatures above that of the surrounding sand. For example, in green (Chelonia mydas) and loggerhead turtle (Caretta caretta) nests, a temperature differential between 0.07°C and 2.86°C has been attributed to metabolic heating [28],[29]. A divergence in measured and predicted temperatures towards the end of incubation would not affect our sex ratio estimate, as metabolic heating should be most pronounced after the TSP. For example, in a direct comparison of sand and nest temperatures at identical depths in a study of N. depressus in the Pilbara, Box (2010) concluded that only two of seven nests showed any detectable sign of metabolic heating [20]. Nevertheless, it would be advisable to include a metabolic heating process (e.g. [30]) into a mechanistic model if traits such as hatchling mortality or incubation duration were of interest.

We did not examine the influence of fine-scale variation in nest orientation, shade or sand properties on nest temperatures in this study, as the Cape Domett rookery is fairly homogenous in these respects. The location of the nest within the beach profile would likely produce the most sand temperature variation due to variable sand moisture and potential shading from vegetation. These parameters can all be modified and mapped spatially using the NicheMapR system (e.g. [8]). Rapid changes in beach topography can make fine-scale modelling of nest shade and exposure challenging, but high-resolution, high frequency satellite and LIDAR imaging could ultimately allow the dynamic rendering of beach temperatures and sex ratios.

3.2 Delineation of the thermosensitive period (TSP)

The precise delineation of the TSP is important for determining the sequence of incubation temperatures that contribute to sex determination [31]. The TSP is the least researched of all TSD parameters and most reptile studies delineate the TSP no more precisely than the middle half to the middle third of development [2],[32]. We showed that the TSP of N. depressus spans a period of about 23% of embryonic development (43%–66%), supporting our hypothesis that the TSP would be less than 33%–50% of the total development period. Hewavisenthi and Parmenter [19] provide the only other estimate of the TSP in N. depressus, concluding that it occurs between days 32 and 40 of incubation at 26°C. Based on incubation temperatures reported in the Hewavisenthi and Parmenter study, and the development rate model developed here (Figure 2), the proportion of development that would have occurred on days 32 and 40 would have been ~44% and ~54%, respectively. Both studies therefore infer that the TSP starts at about the same developmental stage. However, this study infers that the TSP is longer; encompassing approximately 23% rather than 10% of embryonic development shown in the Hewavisenthi and Parmenter study.

All hatchlings (except one) produced from female–male–female (FMF) treatments were female, indicating that either a masculinising switch does not fall within the TSP we identified, or else that our switch periods were too short or not sufficiently cool to trigger the development of males. The switch periods were designed to encompass the same proportions of development for switches in both directions. Instead, development of N. depressus embryos from Cape Domett occurred at slightly faster rates than for embryos collected from the Pilbara [20]. This caused a greater proportion of development to be completed by the start of the FMF switches than for the corresponding male–female–male (MFM) switches.

Switches to 32.5°C over as little as only 9% of development produced a feminising effect (Figure 3; treatment 5). Experiments on the pond turtle (Emys orbicularis) show that bipotential gonads develop into testes in the absence of a feminising temperature trigger during the TSP [33]. In contrast, when this trigger occurs, the enzyme aromatase acts on gonads to stimulate the production of oestrogen and development of ovaries [1],[33],[34]. In N. depressus, we hypothesise that the TSP is similar for triggering testicular or ovarian development, but that the proportion of development at feminising temperatures may be 9% (or less) to trigger ovarian development. Hence, if gonads are to develop into testes, the majority of the TSP needs to occur at a masculinising temperature. This is a possible reason for the lack of masculinisation seen in the FMF back-switch experiments (treatments 7–11), as experimental switch periods may have been too short for the majority of the TSP to occur at the masculinising temperature. Further experimentation is required to test this hypothesis, and its resolution will have important implications for estimating hatchling sex ratios in natural nests that fluctuate between male- and female-producing temperatures.

3.3 TSD parameters and embryonic development rate

Identifying the Tpiv, TRT and embryonic development rate were not major foci of this study, but the ability of our physiological model (Figure 1) to provide good estimates of hatchling sex ratios at Cape Domett (Table 3) suggests that these parameters and the boundaries of the TSP were reasonably accurate. The Tpiv estimated in this study of 29.4°C for the Cape Domett population is similar to the 29.3°C estimated by Limpus (1995, cited in [35]) and the 29.5°C estimated for Eastern Australian populations [18]. It is, however, lower than the 30.1°C estimated for a population of N. depressus from the Pilbara region of Western Australia [20]. Differences in the Tpiv estimates for the two Western Australian populations are concordant with the genetic differentiation reported between summer nesting populations from the Pilbara and winter nesting populations in Northern Australia [36].

We used sex ratio data from field nests to improve our estimates of the Tpiv and TRT, but field data is not an ideal substitute for carefully controlled laboratory experiments. In addition, we have not assessed sex ratios close to the estimated Tpiv (29.4°C). Consequently, constant temperature incubation experiments at temperatures at and near the Tpiv are needed to further refine the TSD parameters for the Cape Domett population. Incubation experiments at higher temperatures would also be highly relevant when considering impacts of global warming, as development rate increases with temperature only up to a certain point, after which it rapidly decreases due to enzyme deactivation [2]. High nest temperatures are a documented cause of mortality in marine turtle nests [12],[37] and well-resolved thermal response curves will be essential for predicting embryonic survival under climate change.

3.4 Implications for population viability and management

Marine turtles have adapted to past climate change and may also adapt to accelerated rates of climate change associated with anthropogenic global warming [38]. Adaptations to counter increasing temperatures include adjusting the timing of nesting, altering nest characters such as shade or depth and altering nesting latitude [13],[39]. As peak nesting at Cape Domett coincides with the coolest sand temperatures, it constrains opportunities for females to temporally adjust peak nesting times to counteract increases in incubation temperatures. Another behavioural characteristic—nest site choice with respect to shade cover—has shown some degree of plasticity in other reptiles, where females select sites of differing shade availability in different climates [40],[41]. The effectiveness of this response may be limited at Cape Domett, as shaded nest sites are scarce. Alternatively, females could dig deeper nests, where sand temperatures are usually cooler and more stable [14], but nest depth is limited primarily by limb length in turtles, meaning that nest depths are relatively fixed relative to those of other reptiles [40]. Similarly, cooler nest locations low on the beach risk flooding and suffocation of embryos.

The heritable variation in TSD traits presents an alternative avenue for adaptation to increasing temperatures. Microevolution of the Tpiv and/or TRT in response to climate change is possible [42],[43], and, unlike behavioural adaptations which alter temperatures experienced by embryos, changes in the Tpiv or TRT would alter the sex ratio produced at a given temperature. In general, data on the heritability of physiological traits in reptiles is scarce and as such provides an important avenue for future research [43]. Mechanistic models such as those developed here readily allow the effectiveness of microevolutionary change to be examined.

The overwhelming consensus that anthropogenic climate change will increase global temperatures means that female biases in primary sex ratios could become more prevalent in marine turtles [44],[45]. Consequently, population viability may decline and in extreme cases, localised extinctions could occur [10]. Projected global warming in Queensland suggests that within 50 years, ratios of one male to four females are possible for many marine turtle rookeries in the Great Barrier Reef [38], and such biases may not be sustainable when combined with other anthropogenic stressors such as habitat modification and loss of rookeries [39],[42]. However, polyandrous mating systems may allow marine turtles to persist under female-biased operational sex ratios (e.g. [46]). Understanding responses to operational sex ratio biases and accurately predicting primary sex ratios will both be important when assessing population viability under a warming climate.

5.1 Ethics Statement

All procedures were reviewed and approved by the University of Western Australia animal ethics committee (UWA animal ethics permit RA/3/100/1145) and the Department of Parks and Wildlife (scientific licence SF008844).

5.2 Study site and egg collection

The N. depressus rookery at Cape Domett occurs on a gently sloping, 1.9 km-long, north-west facing beach in the East Kimberley region of Western Australia, Australia (14.798° S, 128.415° E). Freshly oviposited eggs were collected for transport to a laboratory in Perth, Western Australia on September 20 2012. Approximately, whole clutches (52 eggs) were collected from each of six nesting females (total 296 eggs) and stored in damp vermiculite within a portable refrigerator (model: Engel MT45FP) set to 8°C to suspend embryonic development thereby facilitating egg viability during transport (Harry & Limpus 1989). Eggs were transported by helicopter to Kununurra where they awaited further transport via a commercial airline to Perth. Inadvertently, eggs were held over in Kununurra due to a lack of cargo space, and consequently eggs reached the University of Western Australia about 100 h after oviposition. This time period fell outside the 72-h window of viability defined by Harry and Limpus [47].

In situ incubation temperatures were recorded in 11 nests on the nights of September 16 and 17 in 2012. Four Thermochron® iButtons (Maxim Integrated Products; DS1921H; accuracy ±1°C, precision 0.125°C) were placed at intervals into each nest chamber during oviposition, with the goal being to record temperatures at different depths within the nest chamber. The nests were situated at various locations within the ‘pink sand’ sections of the beach [21] and were marked with a GPS device (Garmin eTrex Vista HCx) and photographed for ease of relocation. A weather station (WeatherHawk, Signature Series 232) was erected midway along the beach on a small dune above the high tide mark. Weather data (air temperature, humidity, wind speed, solar radiation) were recorded at 30-min intervals.

The rookery was revisited 40 days later (25–26 October 2012) and all nests were relocated and excavated. Between 12 and 16 eggs that were closest to each of the four iButtons were removed from each nest (total of 138 from all nests). The depths of iButtons below the sand surface were measured to the nearest 10 mm before removal. All remaining eggs were reburied. The embryos inside the sampled eggs were euthanised by chilling for approximately 24 h. Chilling was achieved on the beach using ice, and eggs were relocated to a portable refrigerator set to 2–4°C within 2–3 h of collection. Chilling was maintained during transport by boat and road back to Kununurra. There, dead embryos were removed from their shells and fixed in 10% buffered neutral formalin for road transport to Perth. Ten of the eleven nests (nest IDs 1, 2, 3, 4, 5, 6, 10, 12, 13 and 15) contained viable embryos; eggs from the remaining nest showed no signs of development.

5.3 Incubation experiments

The six live clutches returned to Perth were either used in ‘back-switch’ experiments designed to delineate the TSP, or else were incubated at a constant temperature close to the estimated pivotal temperature. Ten eggs (2 each from 5 clutches) that were not used in back-switch experiments were incubated at a constant temperature of 30.1°C (Model i180), which was the estimated pivotal temperature for a population of N. depressus from the Pilbara region of Western Australia [20].

For the back-switch experiments, one egg from each clutch was buried in damp sand within one of 48 partially sealed containers. Containers were placed in one of four Steridium incubators (models: i180 and i500) set to either a male-producing temperature (28°C; 24 containers) or a female-producing temperature (32.5°C; 24 containers). During development, each container of 5–6 eggs was switched to the opposite temperature treatment for either 10% or 15% of embryonic development (96–240 h, depending on the time and temperature), after which time they were switched back to their original temperature (Table 1). Each back-switch treatment was replicated within a pair of incubators (one pair of 180 L incubators and one pair of 500 L incubators), and iButtons were placed next to an egg when each container underwent a switch to record the temperatures experienced by the embryos.

We used a development rate function developed for N. depressus by Box (2010) to calculate the day that embryos should have completed 35% of development to hatching when incubated at either 28°C or 32.5°C. Switches then began on this day (or when embryos were further developed—refer Table 1) and concluded when a further 10% or 15% of development to hatching should have been completed at the switch temperature (28°C or 32.5°C). Procedural controls were also included for the earliest switches, in which case embryos incubated at 28°C were switched to a different incubator at 28°C for the appropriate period, and then back again to the original incubator. A procedural control was also applied to eggs incubated at 32.5°C.

Eggs in both the back-switch and Tpiv experiments were weighed at the start of incubation and thereafter every seven days to assess embryonic viability. Unviable eggs (those that lost weight and became discoloured) were removed from incubation boxes. On each occasion, containers were repositioned randomly in the incubators to minimise any impact of subtle temperature gradients within incubators. Containers were checked daily after 35 days of incubation (one week before the earliest expected hatch date) and incubation time was recorded as the time from the start of incubation to when pipping commenced. The four days between oviposition and the start of incubation were not included in the incubation time, as eggs were cooled during this period and embryonic development was suspended [47]. Once pipping commenced, hatchlings were weighed and euthanised by an intracoelomic injection of 0.4 mL sodium pentabarbitone at a concentration of 160 mg/kg. Hatchlings were labelled and stored in 10% buffered neutral formalin for later dissection.

5.4 Identification of hatchling sex

The gonads of marine turtle hatchlings are small (<500 μm) and attached to the kidney [48]. Consequently, the left kidney of each preserved individual was removed. Kidneys were then prepared as paraffin-embedded sections and stained with haematoxylin and eosin for light microscopy [31]. The sex of each specimen was identified based on the differentiation of gonadal medulla and cortex, where males have seminiferous tubules in the medulla and a regressed cortex, and females have a disorganised medulla and a thick, well-developed cortex [49],[50]. Each specimen was classified as male, female or unknown on three separate occasions, without reference to previous assessments, and a repeatability analysis was conducted to determine the reliability of the classification. Any specimens where a gonad was not visible or was unable to be classified were resectioned and reexamined until each individual could be classified as male or female.

5.5 Calculation of a non-linear function for embryonic development rate

The parameters fitted by DEVARA describe the maximum development rate (b₁), its corresponding temperature (b₃) and the temperature at which development rate approaches zero (b₂). Parameters b₄ and b₅ control the asymmetry and steepness of the curve and were fixed at 6 and 0.4, respectively, as is recommended when development rates at extreme incubation temperatures are unknown [51].

5.6 Delineation of the thermosensitive period (TSP)

Using the parameters estimated for the non-linear development rate function defined above, hourly temperature records from the back-switch experiments were converted into developmental increments. These increments were integrated to determine the cumulative proportion of development that was completed on each day. The length of each switch window was thereafter described as a developmental proportion (e.g. 12% of development). If any back-switch regime produced a sex ratio different to that expected from the dominant incubation temperature, then the TSP was assumed to fall within the portion of development where the temperature was switched.

5.7 Determination of the pivotal temperature (Tpiv) and transitional range of temperatures (TRT)

The relationship between incubation temperature and sex ratio was calculated using TSD software version 4.0.3 [http://max2.ese.u-psud.fr/epc/conservation/TSD/index.html] developed by [52]. This program is widely used to compare the fit of up to five threshold models using maximum likelihood (Richards/a-logistic, Weibull, Weibull*, Hill and Hill*) where * indicates the a-logistic or ‘asymmetrical’ version of the model. Data from constant temperature incubation were intended to be used for this analysis, but the program requires that three or more temperatures produce mixed sexes in order for the fit of different models to be compared [52]. Of the five constant incubation temperature regimes applied in this study (Tpiv experiment, male back-switch controls and female back-switch controls), only one produced mixed sexes. Hence, sex ratio and temperature data from field nests were added to the laboratory data to produce a more rigorous estimate of a sex ratio function.

As the TSD software is designed for constant temperature data, temperature records from field nests were converted to constant temperature equivalents (CTEs; [53]). Effectively, the CTE is the temperature above which half of development occurs (i.e. a developmental median) and thereby the relationship between development rate and temperature must first be established. Furthermore, when incubation temperatures are variable, the CTE during the TSP (and not the entire incubation period) is most relevant for fitting the TSD function [31]. Hence, we used our development rate function (Figure 2) to convert the nest temperatures that fell within the TSP into hourly developmental increments. These increments were ranked and integrated to determine the temperature above which half of development occurs (the CTE).

where sr(t) is the sex ratio at a given temperature t, P is the pivotal temperature, S describes the steepness of the transition from male to female producing temperatures and K is a parameter that describes the asymmetrical shape of the function. The fitted values of P, S and K and the corresponding TRT were noted and used as reference points for a further ten models, which were fitted from random subsets of five field nests (Table 2). By using only five of the ten field data points available, we were able to meet our goal of retaining an independent data set to test our mechanistic framework (Figure 1). For each of the ten subsets, there was more than one equation that produced an equally good fit (∆AIC <2); however, the Hill equation consistently produced the lowest AIC value across all subsets. Hence, we used the Hill equation to make a standardised comparison of the TSD equation generated by each of the ten datasets, and used the mean residuals to select a consensus model that could be used to estimate the Tpiv and TRT of the Cape Domett population.

5.8 Sand temperature reconstruction and projection

NicheMapR is an R version of the mechanistic modelling program NicheMapper^TM; [54] and was used to predict hourly sand temperatures by simultaneously solving heat and mass balance equations based on physical properties of beach sand and on regional climate data (e.g. [11],[22]). The parameters of the microclimate model were not ‘tuned’ to the observed sand temperatures. Rather, the sand properties for the Cape Domett rookery were estimated from relevant literature (e.g. [55],[56]) or were measured empirically (Table 5). Percent sand moisture was measured by determining the wet and dry weights of sand samples collected from a range of depths from each field nest during excavation in October. Solar reflectance of sand sampled from two nests in the wavelength range 300–2,100 nm was measured using two spectrometers (Ocean Optics USB2000 for the UV-visible range and NIRQuest for the NIR range) and two light sources (Ocean Optics PX-2 pulsed xenon light for the UV-visible range and HL-2000 tungsten halogen light for the visible-NIR range), all connected with a quadrifurcated fibre optic. The probe on the end of the fibre optic was held within an ocean optics RPH-1 probe holder at a constant angle (45°) and distance from the surface, and each measurement was expressed relative to a Spectralon 99% white reflectance standard (Labsphere, Inc., North Sutton, NH, USA), and weighted by solar irradiance. There was no significant difference in the solar reflectance of the two samples (unpaired t-test; p = 0.959).

Daily maximum and minimum temperatures, relative humidity and wind speeds were collated from data generated by the weather station we installed at the rookery between September 18 and October 24, 2012. Rainfall data during this period were obtained from the nearest Bureau of Meteorology weather station at Wyndham (http://www.bom.gov.au). These weather values were used as inputs into NicheMapR, and sand temperatures were predicted for user-specified depths (range 36–80 cm) using the parameters in Table 2 and assuming 0% shade. Predicted temperatures at a particular depth were then compared to actual nest temperatures measured at the same depth, using the r² statistic.

Historical (1990–2013) climate data for Cape Domett (daily maximum and minimum temperatures, relative humidity, rainfall and solar radiation) were obtained from the Australian Water Availability Project (AWAP). Gridded long-term average wind speed data were obtained from Australian National University Climate software package (ANUCLIM; [57]). Both the AWAP and ANUCLIM databases are derived from interpolated data records from weather stations across Australia. These historical climate data were interpolated at a point approximately 25 km south of the Cape Domett rookery (15.003° S, 128.383° E).

The AWAP climate data, combined with the shade and sand parameters defined earlier, were used as inputs into the microclimate model within NicheMapR to estimate sand temperatures at 50 cm depth for the 23 years from 1990 to 2013. To investigate the influence of climate change on sand temperatures, the 2007 air temperatures from the AWAP database were adjusted in accordance with the Commonwealth Scientific and Industrial Research Organisation’s (CSIRO) projections of future air temperatures for Australia [58]. Under a low emissions scenario, AWAP air temperatures were increased by 0.6°C or 1.8°C for the years 2030 and 2070, respectively. Under a high emissions scenario, air temperatures were increased by 1.5°C or 3.4°C for the same years.

5.9 Model validation

The DEVARA function was used to convert the hourly sand temperatures predicted by NicheMapR (or the actual temperatures measured in nests) into developmental increments for the five nests not used to fit the TSD function. Development was assumed to have started on the date of oviposition; hence, the first developmental increment we calculated was for the hour immediately following oviposition. All the hourly development increments calculated for a particular temperature record were integrated to determine the proportion of development completed on each day and thereby identify the dates that formed the boundaries of the TSP. We calculated a CTE for all temperatures that fell within the TSP (as described earlier) and CTEs were inputted into the TSD function (Figure 4B) to estimate the sex ratio of the embryos sampled at the various depths in each nest. Sex ratios were classified as female (<5% males), mixed (5–95% males) or male (>95% males). Predicted sex ratios were then compared to the actual sex ratios of the sampled embryos, which allowed us to assess how effectively our physiological model estimated sex ratios from hourly temperature records (Figure 1).