The effects of formalin fixation and fluid storage on stable isotopes in rodent hair

Abstract

Stable isotopes are increasingly being used to unlock the wealth of information contained in specimens preserved in museum collections. However, preservation methods that employ formalin may confound ecological interpretations. To quantify the effects of formalin fixation and subsequent fluid storage in ethanol on the isotopic signatures of small mammal hair, we analyzed δ ¹³C and δ ¹⁵N values from specimens of seven rodent species that were sampled repeatedly both before and after varying lengths of formalin fixation (1–11 days) and ethanol storage (1–6 years). We supplemented these data with a 2-week fixation experiment using deer mice (Peromyscus maniculatus) in which no ethanol storage was employed. As expected, preservation in formalin and ethanol had no discernable effect on δ ¹⁵N values. In contrast, specimen δ ¹³C values decreased in a saturating fashion during formalin fixation and over subsequent years of fluid storage in ethanol. On the basis of models that we fit to these time series, we estimate the long-term effect of fixation and storage on δ ¹³C values to be −0.92‰ after 4 years. This biologically relevant shift in δ ¹³C values should be accounted for when inferring the diets of species from fluid-stored museum collections and when comparing across specimens with different preservation histories.

diet reconstruction, ethanol, formaldehyde, isotope ecology, museum collection, small mammals

Natural history collections represent vast and underused repositories of ecological data, including specimens conducive to stable isotope analysis. In the last decade, recognition of the utility and increased accessibility of stable isotope biogeochemistry for studying a species’ niche has opened many new research directions (Newsome et al. 2007; Koch et al. 2009). For mammals, the stable isotope compositions of muscle, blood, bone, and hair now are routinely used to determine the relative trophic position, dietary niche breadth, resource use, and migration patterns of both modern and historically collected specimens (Bearhop et al. 2004; Phillips et al. 2005, 2014; Newsome et al. 2007; Podlesak et al. 2008; Reid et al. 2013). This is because δ ¹³C and δ ¹⁵N values in these tissues reflect aspects of an animal’s diet, with δ ¹³C providing information about the plant type consumed (e.g., C₃ vs. C₄ vs. CAM), δ ¹⁵N capturing information on relative trophic position (Newsome et al. 2007; Koch et al. 2009), and both reflecting aspects of habitat, such as aridity (Newsome et al. 2007). However, the effects of museum preservation techniques on the isotopic signal recorded in small mammal hair, the tissue most readily available in natural history collections, currently are unknown.

Preparation and storage of small mammal specimens vary by museum, collector, and era in which specimens were collected, with methods including, but not limited to, skeletonization, study skin preparation, and fluid preservation (Quay 1974; Simmons and Voss 2009). Fluid preservation in ethanol preserves all body tissues and gross morphology, and most often follows the subcutaneous and intraperitoneal injection of formalin followed by short-term (i.e., days to weeks) storage of specimens in formalin prior to accessioning. Formalin is produced by saturating water with formaldehyde gas to a concentration of 37% by weight and, for specimen fixation, is diluted to 10% as a solution of one part formalin and nine parts water by volume to create an approximately 3.7% formaldehyde solution (Quay 1974; Simmons and Voss 2009). After formalin fixation, specimens then are permanently stored in 70% ethanol (Quay 1974; Simmons and Voss 2009).

Fluid preservation and storage present a unique challenge for stable isotope analysis because formaldehyde methylates proteins with isotopically light carbons (¹²C). Formaldehyde is derived from fossil fuel constituents whose carbons originated from ancient C₃ plants, and thus formaldehyde is relatively ¹³C depleted, ranging from −37‰ to −53‰ (Edwards et al. 2002; Dawson and Siegwolf 2007; Alelú-Paz et al. 2008; Reuss 2012). The post-mortem addition of carbons from a ¹³C-depleted source pool to the proteins in animal tissues manifests as a bulk shift of ¹³C in the tissue towards lower δ ¹³C values (analogous to the Seuss effect; Keeling 1979; Bacastow et al. 1996). In contrast, δ ¹⁵N values should remain unchanged because formalin exposure does not introduce new nitrogen atoms to the tissues.

While ethanol storage has been shown to affect the isotopic signals of blood and muscle tissue via the alteration of lipid content or the retention of ethanol in these tissues, such effects are not expected for bone or hair samples. Specifically, ethanol can dissolve in lipids present in blood and muscle and not be completely removed by drying. For example, Javornik et al. (2019) demonstrated that ethanol storage increased δ ¹³C values in muscle and liver tissue of bears by 0.4‰ and 0.6‰, respectively. This is similar to offsets obtained after the implementation of procedures specifically designed for lipid extraction (e.g., +0.4‰ and +0.8‰ in δ ¹³C for muscle and liver, respectively; Javornik et al. 2019). Hair and bone collagen are unlikely to experience this same ethanol-lipid interaction because they are proteinaceous tissues.

Even though an effect of ethanol on hair is not expected, it cannot be assumed that δ ¹³C values in hair remain constant once a formalin-preserved specimen has been transferred to permanent storage in ethanol. This is because specimens exposed to formalin for several days to weeks prior to ethanol storage can accumulate unreacted formalin in the body cavity that persists after transfer to permanent ethanol storage (Simmons and Voss 2009; Buesa and Peshkov 2012). Even when specimens are rinsed of formalin before ethanol storage, residual external formalin and formalin injected into the body cavity and tissues will remain. Once a specimen is in ethanol, internal formalin can diffuse out of the body cavity following its concentration gradient to become present on the exterior of the specimen (Simmons and Voss 2009). Medewar (1941), working with mammalian soft tissues, suggested that the rate of formalin diffusion through a body is approximately 5 mm per 25 h. Further methylation of proteins by residual formalin thus will proceed slowly over time while specimens are stored in ethanol either until all the residual formalin undergoes reaction or all available reaction sites within the tissues are saturated (Baxter and Walton 1970).

Here we quantified the impact of formalin preservation and subsequent ethanol storage on the isotopic signal recorded by mammal hair. While many studies have investigated the effects of formalin fixation and ethanol storage on the tissues of terrestrial and marine vertebrates and invertebrates, few have leveraged repeated or paired sampling of the same specimens (Edwards et al. 2002; Sarakinos et al. 2002; Baugh et al. 2004). Baugh et al. (2004) is the only study of which we are aware to have considered mammalian hair, and did not pair pre- and post-fixation sampling of the same specimens. In our study, we measured the δ ¹³C and δ ¹⁵N signals of multiple small mammal species whose specimens we sampled repeatedly both before and after varying lengths of formalin fixation (1–11 days) and ethanol storage (1–6 years). We supplemented these data with a 2-week fixation experiment using deer mice (Peromyscus maniculatus) for which no ethanol storage occurred. As expected, we did not detect an effect of formalin fixation or subsequent storage on δ ¹⁵N values. We did, however, observe a strong effect of formalin fixation and storage on δ ¹³C values, the time-course and asymptotic effect of which we estimated using a mixed-effect statistical framework. Our experiments and models provide a means to infer the pre-preservation δ ¹³C values of specimens whose time spent in formalin fixation and ethanol storage is known. Applying such time-dependent correction factors to isotopic data is important when comparing diets across species, evaluating change in diet within a species over time, and when comparing across specimens and/or species that differ in their collection and preservation histories.

Materials and Methods

Specimen collection and hair sampling

Our first dataset of “field-processed” individuals totaled 33 specimens representing 7 species of rodents [number of individuals]: Peromyscus maniculatus [6], Dipodomys merriami [4], D. ordii [4], Microtus longicaudus [4], Onychomys leucogaster [2], Reithrodontomys megalotis [2], and Perognathus mollipilosus [11] (museum catalogue numbers in  Appendix). These were collected with Museum Special snap traps by E. A. Rickart and R. J. Rowe during the summers of 2011, 2012, and 2013, in the Toiyabe, Pine Forest, and Ruby Mountain ranges of Nevada, respectively, as part of an ongoing resurvey of the small mammal communities in the region (Rowe et al. 2010; Rowe et al. 2011; Rowe and Terry 2014). Hairs representing pre-preservation samples were clipped from the sacral region of the lower back between the anterior margin of the innominate and the tail of each specimen shortly after capture. Specimens then were injected with, and submerged in formalin solution (MACRON Chemicals Formalin Solution mixed with water to reach ~10% dilution; VWR International Holdings, Inc., Wilmington, Delaware) on the day of their capture. Specimens collected in 2011 spent a minimum of 1–6 days in formalin, specimens collected in 2012 remained in formalin for a minimum of 2–4 days, and specimens collected in 2013 spent a minimum of 9–11 days in formalin. These submersion times were calculated by subtracting the verbatim capture date from the date of the end of the field season, as recorded in field notes. At the Natural History Museum of Utah, specimens subsequently were sorted by species, rinsed in water, and transferred to long-term storage in 70% ethanol after a varying but unknown number of days. We repeated hair sampling in 2014, 2015, and 2017, from adjacent locations on each specimen’s sacral region, for each of the original 33 specimens.

Our “lab-processed” dataset entailed a controlled fixation study using 10 Peromyscus maniculatus specimens collected by E. A. Rickart in 2017 at the northern boundary of Fish Springs National Wildlife Refuge, Utah (Museum catalogue numbers in  Appendix). For these specimens, pre-preservation hair samples were taken before each specimen was rinsed with tap water and its peritoneal cavity injected with a Formalin solution diluted to 10% with tap water. All specimens then were placed into a glass jar containing 2 liters of the formalin solution. Pre- and post-preservation samples were taken from adjacent locations along the hind-quarters (rump) of each specimen, with post-preservation hair samples taken at 12 h, 24 h, and 14 days of exposure, to capture the range of typical exposure durations observed for field-processed specimens.

Stable isotope analysis

All samples were prepared for stable isotope analysis following methods outlined in Reid et al. (2013). Hair samples were washed three times with water and petroleum ether under sonication, then rinsed twice with water to remove unreacted formalin, ethanol, lipids secreted from the animal, and soil from the environment. Samples then were dried overnight at 60°C and ca. 500 μg of each sample was weighed into individual tin capsules for stable isotope mass spectrometry.

Carbon and nitrogen stable isotope ratios were estimated at the University of California Santa Cruz Stable Isotope Laboratory. Samples were flash combusted using a Carlo Erba 1108 elemental analyzer (Carlo Erba, Milan, Italy) converting them to CO₂, N₂, and H₂O. The gaseous sample then was passed through a desiccant column to absorb H₂O, leaving CO₂, and N₂, to pass into a ThermoFinnigan Delta Plus XP isotope ratio mass spectrometer (Thermo Fisher Scientific Inc., Waltham, MA, USA). Values were corrected for mass and drift using calibrated in-house standards (gelatin and acetanilide). The mean standard deviations for replicates of these in-house standards across all runs ranged from 0.03‰ to 0.12‰ for δ ¹³C and 0.04‰ to 0.08‰ for δ ¹⁵N values. Corrected isotopic values are expressed relative to the international standards PDB (Pee Dee Belemnite) for δ ¹³C and Air for δ ¹⁵N. Isotope ratios are reported in parts per thousand deviation from the standards by δ (‰) = [(R_sample/R_standard) − 1] × 1,000, where R is the ratio of the heavy isotope to the light isotope. Due to a machine failure that compromised the δ ¹⁵N values for some “lab-processed” samples, these samples were re-run when remaining sample amounts allowed. This resulted in a subset of samples having two uncompromised values of δ ¹³C and δ ¹⁵N. For this subset of specimens, we used the mean of each specimen’s run 1 and run 2 values for its value of δ ¹³C and used run 2 for its value of δ ¹⁵N. (The mean pairwise difference between duplicates was 0.001‰ for δ ¹³C and 0.10‰ for δ ¹⁵N.) The lack of remaining sample resulted in the loss of 9 (of 153) values for δ ¹⁵N.

Statistical analyses

We used a nonlinear mixed-effects framework to fit and estimate the parameters of three different models, each describing the time-course of δ ¹⁵N or δ ¹³C values resulting from formalin fixation and subsequent ethanol storage. Our first model represented the null hypothesis that isotope values (either δ ¹⁵N or δ ¹³C) are not affected by fixation or subsequent storage in ethanol. This model (subsequently referred to as the null model) is written for δ ¹³C as

$${\displaystyle \begin{array}{l}{\displaystyle {\delta }^{13}C(t)=1}\end{array}}$$

(1)

and describes a specimen’s isotopic value starting and remaining at a constant value I [i.e., I = δ ¹³C(t = 0)], with any observed variation over time reflecting some combination of sampling and observation error. We expected this model to perform well at describing the time-course of δ ¹⁵N values and to perform poorly at describing the time-course of δ ¹³C values.

The basis of our two alternative models was the monomolecular function, f(t) = A(1 – w e^-kt), where f(t) represents the concentration of a reactant at time t, A represents the reactant’s asymptotic concentration, k represents the reaction rate at which the asymptote is approached, and w represents the net amount of change per unit of the reactant’s asymptotic concentration (Weber 1891; Brody 1945; Koya and Goshu 2013). This function commonly is used to describe the rate of irreversible chemical reactions and is based on the principal of diminishing increments (France et al. 1996).

The first of our alternative models (subsequently referred to as the biphasic model) describes the time-course of a specimen’s isotope values as being dictated by a two-step process of differing fixation and subsequent storage effects. For δ ¹³C, this model is written as

$${\displaystyle \begin{array}{l}{\displaystyle {\delta }^{13}C(t_f,t_s)=I+A(1-e^{k_f{t}_{f+}{k}_s{t}_s})}\end{array}}$$

(2)

where I represents the specimen’s initial δ ¹³C value [i.e., I = δ ¹³C(t_f = 0, t_s = 0)], k_f represents the rate of change during the t_f number of days a specimen spends in formalin, k_s represents the rate of change during the subsequent t_s number of days a specimen is stored in ethanol, and A represents the asymptotic magnitude of change that a specimen’s value can exhibit (dictated by the specimen’s potential methylation sites). In the context of the monomolecular function, we assumed w = 1 for both fixation and storage, meaning that all potential methylation reactions indeed were achievable. The biphasic model therefore reflects a process by which a constant proportion of the remaining un-methylated proteins are methylated over time, saturating to a value of I + A at a rate of k_s per day in storage after the initial effect of fixation at a rate of k_f per day has been incurred.

The second of our alternative models (subsequently referred to as the monophasic model) treated fixation and storage as representing a single process of continual formalin exposure and did not distinguish between time spent in fixation versus subsequent storage. Representing the total time of preservation as $t_p=t_f+t_s$⁠, the model describes the time-course of isotope values by

$${\displaystyle \begin{array}{l}{\displaystyle {\delta }^{13}C(t_p)=I+A(1-e^{k_p{k}_p})}\end{array}}$$

(3)

(written for δ ¹³C), where I and A retain the same interpretation as before and k_p reflects the rate at which values change after preservation is first initiated. Although the biphasic model may make more intuitive sense given the abrupt nature of a specimen’s transition from formalin to ethanol, this monophasic model has one less parameter needing to be estimated from data and may be equally parsimonious in its description of the underlying process if methylation continues to occur at a similar rate during ethanol storage due to residual formalin within or on the surface of the specimen.

In fitting the three models to our δ ¹⁵N and δ ¹³C data, we combined both our field-processed and lab-processed specimen data together into one dataset and accounted for the repeated sampling of the same specimen by treating initial values, I, as a random effect at the level of the individual specimen. That is, we allowed individual specimens to differ in their initial value from an estimated population-level mean initial value. (Our null model therefore represents nothing more than the mean and variance of all samples.) The other parameters (A, k_f, k_s, k_p) of the biphasic and monophasic models, in which an individual specimen’s value was allowed to change over time, were treated as fixed effects.

We compared the relative support for the three models using both the sample-size corrected Akaike Information Criterion, AIC_c (Hurvich and Tsai 1989), and the Bayesian Information Criterion, BIC (Schwarz, 1978). These information criteria penalize a model’s likelihood (i.e., its “goodness-of-fit” to the data) by its parametric complexity (i.e., its number of parameters given the sample size of the data) to permit evaluations of model performance across models of differing complexity while avoiding the choice of over-fit, and hence unnecessarily complex, models (sensu Occam’s razor, Burnham and Anderson 1998). AIC_c is preferred over AIC (Akaike 1973) when sample sizes may be small relative to a model’s complexity and converges on AIC as the sample size increases (Burnham and Anderson 1998). BIC is more “conservative” than AIC, penalizing models of greater complexity more heavily (Burnham and Anderson 1998). The model with the lowest IC is considered the best-performing model, with alternative models having values within two IC units (e.g., ΔAIC_c < 2) considered to have equivalent support. In such cases, the simpler model is preferred (Burnham and Anderson 1998).

Determining a specimen’s unknown initial δ ¹³C value

Although the fitting of the above biphasic and monophasic models is appropriate for the estimation of their parameters from time-series data, neither model is applicable directly for using these parameter estimates to hindcast a museum specimen’s unknown initial pre-preservation value from its observed post-preservation value. More specifically, use of the rearranged version of the biphasic model, for example, only is appropriate if the effect rates and time periods neither of fixation nor storage were large or long enough to have caused a specimen’s value to reach its asymptote (to within analyzer uncertainty; Fig. 1A); the initial value would be overestimated in cases where a specimen’s value reached its asymptote during storage (Fig. 1B) or fixation (Fig. 1C). Instead, a specimen’s initial value may be determined by describing the preservation process in reverse using the mechanistically equivalent contingent model,

Fig. 1.

The conditional nature of hindcasting a specimen’s initial isotopic value given its observed, post-preservation value, here illustrated for the biphasic model (equations 2, 4, and 5). (A) If neither the effect rates (k_f and k_s) nor time periods (t_f and t_s) of fixation or storage were large or long enough to have caused a specimen’s isotopic value to reach its asymptotic value (solid circle) to within assumed precision, then its initial value (open circle) may be estimated by a simple rearrangement of equation 2. However, if the specimen’s isotopic value reached its asymptote during (B) storage in ethanol or (C) during fixation, then a simple rearrangement of equation 2 will lead to an overestimate of the specimen’s initial value. In the latter two cases, and as described in the main text, correctly estimating a specimen’s initial value requires use of equations 4 and 5. Dashed grey-black trajectories beneath the asymptotic value reflect the erroneous reverse application of the forecasting model to the initial value. The equivalent logic applies to the monophasic model (equation 3).

$${\displaystyle \begin{array}{l}{\displaystyle {\delta }^{13}{C}_{initial}={\delta }^{13}{C}_{obs}-A(e^{k_f}{}^{{{\displaystyle \widehat{t}}}_f}{}^{+k_s{{\displaystyle \widehat{t}}}_s}-1)}\end{array}}$$

(4)

where

$${\displaystyle \begin{array}{l}{\displaystyle ({{\displaystyle \widehat{t}}}_f,{{\displaystyle \widehat{t}}}_s)=\left\{\begin{array}{ll}& {\displaystyle (-t_{f,}-t_s)if{t}_f<t_{f(max)}and{t}_s<t_{s(max)}}\\ & {\displaystyle (-t_{f,}-t_{s(max)})if{t}_f<t_{f(max)}and{t}_s\ge t_{s(max)}}\\ & {\displaystyle (-t_{f(max),}0)if{t}_f\ge t_{f(max)}}\end{array}\right\}}\end{array}}$$

(5 a-c)

Under this model, the effect of fixation only is reversed if the time period that a specimen spent in fixation was less than the fixation time period that would have been necessary for it to have reached its asymptotic value. This hypothetical time period, t_f(max), may be determined using estimates of A and k_f as $t_{f(max)}=In(p({\delta }^{13}{C}_{obs})/A)/k_f$⁠, where p(δ ¹³C_obs) refers to the precision to which the specimen’s post-preservation δ ¹³C value was measured. Similarly, the effect of post-fixation storage is only reversed if the total time that a specimen spent in fixation and subsequent storage was less than the hypothetical time period that would have been necessary for the storage effect to have caused a specimen’s value to reach its asymptote, calculated as $t_{s(max)}=In(p({\delta }^{13}{C}_{obs})/(A{e^{k_f}}^{t_f})/k_s$⁠. Note that the calculation of t_p(max) for the monophasic model is equivalent to the calculation of t_f(max) with k_f replaced by k_p.

Finally, although equations 4 and 5 (and their monophasic model equivalents) may be used directly to obtain a point prediction for a specimen’s initial pre-preservation isotope value, knowledge of the uncertainty associated with this hindcast estimate often is of equal importance. We therefore wrote a simulation-based implementation of equations 4 and 5 (and their monophasic model equivalents) to characterize this hindcast uncertainty by using values sampled from the uncertainty distribution of the corresponding statistical model’s data-estimated parameters (Mandel 2013). The method involves the iterative sampling of parameters to make repeated model hindcasts, the percentiles of the resulting distribution of which are used to construct a 95% prediction (rather than confidence) interval for the point prediction. In the context of the biphasic model, we did so by assuming the parameter estimates of A, k_f, and k_s to be multivariate normally distributed as defined by their estimated variance-covariance matrix. The covariance of δ ¹³C_obs with A, k_f, and k_s was set to zero due to their statistical independence, while the variance of δ ¹³C_obs was assumed to be equivalent to the variance of the replicate in-house gelatin controls that were run when the specimen’s hair sample was isotopically analyzed. The equivalent assumptions apply to the monophasic model. We wrote the functions for applying our simulation-based hindcasting models in such a way that they can provide either (i) hindcast predictions of the initial values of an arbitrary number of observed specimens, or (ii) a hindcast time-series extending from the observed to the predicted initial value of a single specimen.

All analyses were carried out in R (v. 3.5.2, R Core Development Team 2018) and made use of the following packages: “lme4” (Bates et al. 2015), “nmle” (Pinheiro et al. 2018), and “MASS” (Venables and Ripley 2002). All data and code, including functions for applying our simulation-based hindcasting models to new data, are available at FigShare (DOI: 10.6084/m9.figshare.12664859) and at https://github.com/marknovak/d13C_hindcast. The data also are provided in Supplementary Data SD1 and SD2.

Results

Stable isotope ratios

Visual inspection of the unmodeled δ ¹⁵N stable isotope ratios for the 33 field-processed and 10 laboratory-processed individuals suggested that their isotope ratios did not, on average, change over time, either during fixation or subsequent storage (Fig. 2). This was confirmed by our statistical analysis in that the model that performed best at describing the time course of δ ¹⁵N values was the null model. This model described a specimen’s true isotopic value as remaining constant over time, with all observed variation occurring due to sampling and observation error and differences among specimens in their initial, pre-fixation value (Table 1). As estimated under the null model, δ ¹⁵N values exhibited a population-level mean of 9.47‰ (± 0.57 SE, t-value = 16.47) and a standard deviation of 3.75‰ (0.47 residual) across specimens. By comparison, the monophasic and biphasic models provided initial δ ¹⁵N estimates of 9.37‰ (± 0.58 SE, t-value = 16.04) and 9.35‰ (± 0.59 SE, t-value = 15.93), respectively, and provided estimates indicating negligible temporal effects. That is, for the monophasic model: A = 0.22‰ (± 0.18 SE, t-value = 1.21) and k_p = 1.53 × 10⁻³ ‰/day (± 3.77 × 10⁻³ SE, t-value = −0.40); and for the biphasic model: A = −0.11‰ (± 0.36 SE, t-value = −0.30), k_f = 0.11‰/day (± 0.19 SE, t-value = 0.59), and k_s = 3.10 × 10⁻⁴‰/day (± 6.63 × 10⁻⁴ SE, t-value = 0.47).

Table 1.

Information-theoretic comparison of the relative performance of our three models in fitting the time-course of δ ¹⁵N and δ ¹³C values.

Isotope	Model	df	ln(ℒ)	Δln(ℒ)	AICc	ΔAICc	AICc weight	BIC	ΔBIC	BIC weight
δ 15N	Null	3	−210.5	0	427.2	1.6	0.25	435.9	0	0.929
	Monophasic	5	−208.7	1.8	427.8	2.2	0.19	442.2	6.3	0.041
	Biphasic	6	−206.5	4.0	425.6	0	0.56	442.8	6.8	0.031
δ 13C	Null	3	−197.3	0	400.7	84.4	<0.001	409.7	78.6	<0.001
	Monophasic	5	−152.9	44.3	316.3	0	0.930	331.0	0	0.930
	Biphasic	6	−153.0	44.3	318.6	2.3	0.070	336.2	5.2	0.070

df refers to its total number of free parameters (including estimates of uncertainty). ln(ℒ) refers to the model’s log-likelihood. The model with the lowest AIC_c or BIC (e.g., ΔBIC = 0) is considered the best-performing model, with alternative models within two IC units (e.g., ΔBIC < 2) considered to have equivalent support. A model’s weight may be interpreted as the conditional probability that it is the best-performing model.

Fig. 2.

The effect of storage in ethanol on A) δ ¹³C and B) δ ¹⁵N ratios of repeated samples from field-processed specimens prior to and after formalin fixation. Gray-scale of points indicate minimum time spent in formalin prior to storage in ethanol. The effect of formalin fixation on (C) δ ¹³C and (D) δ ¹⁵N ratios for repeated samples from laboratory-processed specimens which spent no time in ethanol.

For δ ¹³C, visual inspection of the unmodeled δ ¹³C stable isotope ratios for the 33 field-processed and 10 laboratory-processed individuals suggested that ratios became lower over time (Fig. 2). This was also confirmed by our statistical analysis. The model that performed best at fitting the time course of values was the monophasic model. This model described a specimen’s isotope value as declining with respect to the total time of preservation regardless of the differentiation between fixation and subsequent ethanol storage (Table 1, Fig. 3). That is, there was strong support in our data for an effect of fixation and storage on δ ¹³C values, but little support for this effect being better described by a biphasic process. As estimated using the monophasic model, initial δ ¹³C values I exhibited a population-level mean of −22.64‰ (± 0.04 SE, t-value = −498.99) and a standard deviation of 2.04‰ (0.14 residual) across specimens. δ ¹³C values then were inferred to decline at an estimated rate k_p of −3.13 × 10⁻³ ‰/day (±1.08 × 10⁻³ SE, t-value = −2.90) of preservation before reaching an asymptotic total decline A of −0.92‰ (± 0.07 SE, t-value = −13.11). The biphasic model produced near-equivalent point estimates for I and A, providing a population-level mean I of −22.64‰ (± 0.01 SE, t-value = −2466.38) and a standard deviation of 2.02‰ (0.14 residual) across specimens, an asymptotic total decline A of −0.92‰ (± 0.01 SE, t-value = −101.03), and rates of decline for fixation k_f and ethanol storage k_s of −4.49 × 10^–4 ‰/day (± 7.18 × 10⁻³ SE, t-value = −0.06) and −3.21 × 10⁻³‰/day (± 1.22 × 10⁻³ SE, t-value = −2.62), respectively. Of direct relevance to our method for quantifying hindcast uncertainty, the estimated correlations among the parameters of the monophasic model were

Fig. 3.

The predicted time-course of a specimen’s δ ¹³C values from its initial value (denoted by solid point) to its estimated asymptotic value (denoted by dashed line) as described by the A) monophasic and C) biphasic models. For illustrative purposes, we here assume a 14-day period of fixation followed by 4 years of storage in ethanol, as was typical for our field-processed specimens. In A) and C), upper and lower bounds delineate the 95% prediction interval which reflects both the estimated variances and covariances of the parameter point estimates. The magnitude of uncertainty depicted in B) and D) refers to the difference between the upper and lower bounds of the prediction interval for the monophasic and biphasic model, respectively. The uncertainty in the initial value reflects the uncertainty from both the fixed effect (population-level mean) and the random effect among individuals. Prediction uncertainty then increases due to uncertainty and covariances among parameter estimates, and then declines as the more precisely estimated asymptotic effect of preservation is approached.

I	A
A −0.656	−
kp −0.086	−0.363

while for the parameters of the biphasic model, they were

I	A	kf
A −0.001	−	−
kf −0.001	−0.001	−
ks −0.003	−0.038	−0.079

In Fig. 4A, we used the point and covariance estimates described above to illustrate the hindcasting of a hypothetical specimen’s initial, pre-preservation δ ¹³C value. In so doing, we set this hypothetical specimen’s measured value to the average δ ¹³C value of all field-processed specimens at their final measured time point, which was −23.8‰, and set the uncertainty of this value to the most representative standard deviation observed for the replicate in-house gelatin controls that corresponded to these field-processed specimens, which was 0.04‰. We assumed that the time periods of fixation and storage were 10 days and 10 years, respectively, and an analyzer precision of two decimal places. (We also evaluated the effects of assuming a precision of one or three decimal places.) Fig. 4B provides a comparison of the true initial values of all field-processed specimens to their hindcast initial values as estimated using their observed final values, actual times spent in preservation, and the standard deviation of replicate in-house gelatin controls within the run in which each sample was analyzed.

Fig. 4.

A) The hindcast time-course of a hypothetical specimen’s δ ¹³C values—given an observed, post-preservation value measured after some known time period of preservation—assuming the best-performing monophasic model and its parameters and covariances as estimated from our data. For illustrative purposes, we here assume a 10-day period of fixation followed by 10 years of storage in ethanol. The specimen’s post-preservation observed value is here assumed to be equal to the average final δ ¹³C value of all field-processed specimens (denoted by solid point), with a standard deviation equaling the standard deviation of replicate in-house gelatin controls. Light and dark grey bands, respectively, indicate the 95% and 68% confidence intervals. The transition from lighter to darker gray regions after a preservation time of 1443 days (3.95 years) reflects the point at which the specimen’s value was inferred to have reached its observed (and asymptotic) value to within a precision of two decimal places. B) Comparison of the actual and predicted initial δ ¹³C values (± 95% confidence interval) of the field-processed specimens. Predictions were made for each specimen using its observed final δ ¹³C, its time spent in preservation, and the standard deviation of replicate in-house gelatin controls within the run in which the final sample was analyzed.

Discussion

Stable isotope ratios recorded in a specimen’s tissues offer a tool for unlocking behavioral and ecological data from previously collected specimens, providing proxy measurements both for resource-based and habitat-based niche dimensions (Newsome et al. 2007). We hypothesized that preservation of museum specimens with formalin generates a consistent directional bias in δ ¹³C signals of small mammal hair that could alter these biological inferences, because formalin is derived from fossil fuel constituents that are ¹³C depleted (Keeling 1979; Bacastow et al. 1996). We observed this expected depletion in δ ¹³C values over time using repeated samples of small mammal hair exposed to formalin and stored in ethanol, estimating an asymptotic effect of 0.92 ‰. This asymptotic effect is expected to be incurred after 3.95 years when δ ¹³C values are measured to a precision of two decimal places. A precision of 1 decimal place would result in an estimated time period of 1.95 years, while a precision of 3 decimal places would result in an estimated time period of 5.96 years. As expected, we did not detect an effect of fixation or storage on δ ¹⁵N values.

Because carbon isotopes are well-known to reflect the composition of an animal’s diet with respect to the photosynthetic pathway of plant resources (e.g., C₃ vs. C₄ vs. CAM), a shift of up to −0.92‰ in δ ¹³C values is a biologically relevant change that impacts interpretation both of generalist and specialist small mammal dietary niches, especially in the context of dietary mixing models (Phillips et al. 2005, 2014). Correcting for this formalin-induced bias in δ ¹³C data thus is necessary for robust reconstruction of the isotopic niche of populations and species, and critical when integrating data from specimens representing different museum collection types (e.g., fluids, study skins, and skeletons) or storage durations.

Despite the prevalence of formalin fixation and ethanol storage for museum curation for many groups of vertebrates and invertebrates, studies on the isotopic effects of such techniques are rare for fluid-preserved mammals. What work has been done has focused primarily on non-mammalian systems. This is not entirely surprising because fish and herpetological collections always have been prepared by fluid preservation, while mammal collections include dry preservation techniques of skins and skeletons as well. Several laboratory-based studies of marine and freshwater invertebrates and vertebrates have demonstrated that formalin fixation alters δ ¹³C values in various tissues by ca. −1 to −3‰. Edwards et al. (2002) described a shift in fish muscle tissues fixed in formalin, with a mean shift of −2‰ after 10 days of exposure that remained constant out to 190 days in formalin. Edwards et al. (2002) also evaluated the shift in δ ¹³C values in fish tissues collected from long-term storage in ethanol (12–15 years), observing an average difference of −0.8 ‰ between formalin-fixed and unpreserved specimen samples. Unfortunately, these formalin-fixed and unpreserved samples were not obtained from the same individuals, potentially confounding estimates of a formalin effect with baseline differences among individuals. The significance of the random effects components of our models shows that tracking isotopic shifts within an individual is important. The use of unpaired fixed and unfixed samples will reduce not only statistical power but also can alter explanatory inferences when sample sizes are insufficient to overcome individual variation.

Sarakinos et al. (2002) undertook an experiment similar to that of Edwards et al. (2002) with freshwater invertebrates and fish muscle tissue, assessing the effects both of formalin and ethanol independently. However, they did not describe methods for removing formalin or ethanol from the tissues of specimens. It therefore is unclear whether they measured isotopic change in the specimens themselves because they instead may have measured residual preservative left on the specimens. Surprisingly, both Edwards et al. (2002) and Sarakinos et al. (2002) found that δ ¹⁵N were reduced by ca. −0.5‰, which is not explained by the biochemistry of formalin fixation (Dapson 2007; González-Bergonzoni et al. 2015). Our results using repeated samples uncover no detectable trend in δ ¹⁵N for rodent hair.

With respect to mammals, work by Baugh et al. (2004), which focused on the effect of formalin fixation and ethanol storage on the δ ¹³C and δ ¹⁵N signals in the hair of shrews, concluded that formalin fixation and ethanol storage had a negligible effect. However, preserved and unpreserved samples were not paired, nor was a time-series of samples evaluated, leading to a lack of statistical power and other potentially confounding issues for estimating initial δ ¹³C and δ ¹⁵N values, as discussed above. Paired sampling of the same individuals is essential for two main reasons. First, fixation is a chemical process and may occur at different rates in different individuals due to individual variation in each animal’s tissues, and is a physical process with respect to diffusion and may therefore vary among individuals due to different attributes (e.g., body size). Second, when formalin-fixed and untreated samples are drawn from specimens of different populations, the mean difference between formalin-fixed and untreated samples may be confounded by other effects (e.g., Seuss effect and variation in resource base) which were not controlled between the two population samples, even when the amount of individual variation within each population is tightly constrained.

Other, non-statistical issues also persist. For example, the sample preparation methods used by Baugh et al. (2004) on shrews consisted of rinsing the specimens in 70% ethanol prior to clipping hair, which is insufficient to remove lipids from the hair (Reid et al. 2013). This will have been especially problematic because shrews are known to be relatively oily mammals, secreting oil from their sebaceous glands that then coats their hair (Eadie 1938). As a consequence, the observed differences between preserved and unpreserved samples likely did not capture the actual fractionation effect of formalin fixation and could not account for error derived from differences in the isotopic signal between lipids and hair protein (Tieszen et al. 1983; Lee-Thorp et al. 1989; Boecklen et al. 2011). Furthermore, lipids tend to be depleted in heavy isotopes, such that failing to remove them from hair would confound the isotopic effect of fixation (Boecklen et al. 2011; Rioux et al. 2019). Recent work by Javornik et al. (2019) explored this lipid effect in more detail using paired samples of brown bear (Ursus arctos) liver and muscle stored in ethanol. While not directly comparable to our study because they did not explore the effects of formalin on tissue isotope ratios, they did report shifts in δ ¹³C values in both tissue types following storage in ethanol. However, samples taken from ethanol storage were not lipid extracted prior to isotopic analysis and the storage offset values they recovered are similar to those they obtained when lipid extraction procedures are employed on tissues not stored in ethanol (+0.4‰ and +0.8‰ for muscle and liver post-storage, respectively). Lipid extraction is important because the fractionation of δ ¹³C that occurs during lipid synthesis results in an average of 6–8‰ difference between proteins and lipids, thus the mixing of protein and lipid isotopic signals confounds isotopic interpretations of tissues with differing lipid content (Post et al. 2007). In addition, the positive correlation that exists between δ ¹³C values and the lipid content of a tissue results in a lipid bias, such that tissues with high lipid concentration could have δ ¹³C values 3–4‰ more negative than lipid-extracted tissues (Post et al. 2007). Thus, it is not clear if the offset recovered by Javornik et al. (2019) is due to ethanol storage causing alteration of the carbon isotope signal in biochemically unknown ways, or if storage in ethanol resulted in partial lipid extraction. While the lipid content of mammal hair may reasonably be expected to be lower than that of muscle or liver tissue, there is evidence that lipids present in and on mammal hair are sufficient to influence both δ ¹³C and δ ¹⁵N values (Rioux et al. 2019) The findings of Rioux et al. (2019) emphasize the importance of undertaking lipid extraction on hair, despite the fact that it is a proteinaceous tissue. All samples in our study were lipid extracted to eliminate the potential of mixing signals measured in lipids and proteins, to wash away any unreacted formalin, and to ensure no other external sources of ¹³C (e.g., cellular hair roots) were influencing the observed offset between samples treated and not treated with formalin.

It is important to note that δ ¹³C shifts due to formalin fixation in different tissues (hair, liver, blood, and bone) are not likely to be directly comparable with one another, and any inferences would be even more tenuous if they were applied across broad taxonomic groups. This is because protein composition differs among tissues within individuals, and the 1°, 2°, and 3° structure of proteins may differ within the same tissues across broad taxonomic groups (Marshall et al. 1991; Fraser and Parry 2011; Geiger et al. 2013). These compositional differences result in variable numbers of methylation sites available for formaldehyde molecules to adhere to, which in turn will lead to different degrees of methylation and consequently different asymptotic change in δ ¹³C values for each tissue type (parameter A in our models). Fish, freshwater invertebrate, and mammal tissues therefore require taxon-specific and tissue-specific estimates. That said, differences in the number of methylation sites need not necessarily alter the per site rate of methylation (parameter k in our models). We therefore hypothesize that future estimates of methylation rates in different taxa may be similar to those inferred by our study.

Our own study is associated with caveats as well. First, our method is not applicable to all tissue types or taxa, and because records of the process and duration of preservation and storage are not always clear from the information associated with individual specimens that one can obtain from specimen tags and/or online data repositories, may not be appropriate for rodent hair samples without any storage information. Our approach therefore relies on access to original survey field notes to ascertain the collection date of specimens and thus to estimate the start point of fluid preservation. Field notes often do contain a record of the duration of time specimens spent in formalin, or provide sufficient information to estimate minimum duration in formalin, as done in our study. But the uncertainty associated with the use of minimum duration estimates will increase the uncertainty associated with a specimen’s corrected value. In addition, formalin solutions sometimes are prepared in the field using powdered paraformaldehyde and stream water, as opposed to pre-prepared formalin solutions that often contain methanol and use purified water. Indeed, even throughout the duration of a field season, powdered paraformaldehyde, liquid formalin, or water may be added to containers used to store specimens (Eric A. Rickart, University of Utah Museum of Natural History, personal communication, November 2017). The consequence is that the actual concentration of formalin used in the field, both for injection and storage, may both be unknown and variable over a single collection period, introducing an additional source of uncertainty into stable isotope values which our method of hindcasting cannot account for.

Additional controlled studies on mammals clearly are needed to evaluate the degree to which an isotopic offset depends on the specific volume and concentration of formalin that is injected into a specimen at the time of preservation relative to a specimen’s size. For example, Simmons and Voss (2009) indicated that bats less than 100 g only may require 3–4 days in formalin to become completely fixed, while bats greater than 100 g require at least a week. Their field methods differed substantially from ours, however, because bats’ abdominal walls were cut open for the explicit purpose of allowing formalin to enter the abdominal cavity (Simmons and Voss 2009). For the specimens we analyzed, formalin was injected directly into the leg muscles and body cavities, but body cavities either were left closed or minimally cut (punctured) to obtain liver samples. Unfortunately, the volume of formalin injected into the body seldom is recorded. Potential differences in body size and formalin treatment methods among individuals (or taxa) are not accounted for in our data and statistical models beyond the inclusion of a random effect of individual for the intercepts, and therefore cannot be accounted for in the current formulations of our hindcasting models. Work incorporating such additional individual- and species-level variables likely will offer useful insight and potential for generalization.

Future work should focus on resolving the functional form of the isotopic effects of time spent immersed in formalin prior to ethanol storage, and the degree to which this functional form may be tissue-specific. Key questions include whether the time to saturation of methylation sites is longer or shorter for different tissues. For example, bone collagen may come in contact with formalin more slowly as it diffuses through the body, thus may be less effected as the isotopically light carbon is integrated into the collagen portion of the bone.

Developing tissue-specific correction factors for specimens preserved with formalin is important to extend the use of existing collections of old specimens for new biogeochemical analyses. Our work evaluated the effects of initial fixation and the duration of time over which residual formaldehyde continues to alter the stable isotope signal of mammalian hair samples stored in ethanol. Our results demonstrate that a correction of δ ¹³C values should be applied to hair samples from specimens that have been fixed with formalin and stored in ethanol in order to robustly reconstruct the isotopic niches of individuals, populations, and species. Our results also suggest that no such correction is required for δ ¹⁵N values given specimen treatment and sample preparation methods that are equivalent to those used here. We recommend that ecological analyses using specimens known to have been fixed in formalin and stored long term in ethanol account for this shift in carbon isotope values.

Supplementary Data

Supplementary data are available at Journal of Mammalogy online.

Supplementary Data SD1.—Metadata information for Supplementary Data 2.

Supplementary Data SD2.—Data table of stable isotope ratios for all samples used in the analyses.

Acknowledgments

We thank R. J. Rowe and E. A. Rickart for collecting pre- and post-preservation samples in the field, and providing access to museum specimens for subsequent re-sampling and space and specimens for running our formalin experiment. We thank Bridget Regan and Morah West for assistance with preparation of hair samples for stable isotope analysis, and D. Andreasen and C. Carney at the UC Santa Cruz Stable Isotope Lab. We also thank R. J. Rowe and members of the Terry Lab for many helpful discussions and comments on previous drafts that helped strengthen this project. This research was supported by funding from NSF DEB-1457500 to RCT.

Appendix

Table A1.

Specimens used in this study

Museum	Catalog number	Locality	Latitude	Longitude	Species
UMNH	34032	United States, NV	38.96217	−117.273	Perognathus mollipilosus
UMNH	34275	United States, NV	39.40004	−117.039	Perognathus mollipilosus
UMNH	34276	United States, NV	39.40004	−117.039	Perognathus mollipilosus
UMNH	34277	United States, NV	39.40004	−117.039	Perognathus mollipilosus
UMNH	34278	United States, NV	39.40004	−117.039	Perognathus mollipilosus
UMNH	34355	United States, NV	39.37684	−117.006	Perognathus mollipilosus
UMNH	34357	United States, NV	39.37684	−117.006	Perognathus mollipilosus
UMNH	34359	United States, NV	39.37665	−117.01	Perognathus mollipilosus
UMNH	34360	United States, NV	39.37684	−117.006	Perognathus mollipilosus
UMNH	34361	United States, NV	39.37684	−117.006	Perognathus mollipilosus
UMNH	34362	United States, NV	39.37731	−117.009	Perognathus mollipilosus
UMNH	35360	United States, NV	41.67319	−118.598	Microtus longicaudus
UMNH	35417	United States, NV	41.59412	−118.422	Dipodomys merriami
UMNH	35418	United States, NV	41.59412	−118.422	Dipodomys merriami
UMNH	35419	United States, NV	41.59412	−118.422	Dipodomys merriami
UMNH	35420	United States, NV	41.59459	−118.427	Dipodomys merriami
UMNH	35431	United States, NV	41.59412	−118.422	Dipodomys ordii
UMNH	35432	United States, NV	41.59412	−118.422	Dipodomys ordii
UMNH	35433	United States, NV	41.59348	−118.425	Dipodomys ordii
UMNH	35434	United States, NV	41.59348	−118.425	Dipodomys ordii
UMNH	35437	United States, NV	41.59412	−118.422	Onychomys leucogaster
UMNH	35438	United States, NV	41.59348	−118.425	Onychomys leucogaster
UMNH	35644	United States, NV	41.71322	−118.768	Reithrodontomys megalotis
UMNH	35645	United States, NV	41.71322	−118.768	Reithrodontomys megalotis
UMNH	35723	United States, NV	41.6692	−118.69	Microtus longicaudus
UMNH	35724	United States, NV	41.6692	−118.69	Microtus longicaudus
UMNH	35725	United States, NV	41.6692	−118.69	Microtus longicaudus
UMNH	36542	United States, NV	40.59728	−115.379	Peromyscus maniculatus
UMNH	36543	United States, NV	40.59728	−115.379	Peromyscus maniculatus
UMNH	36544	United States, NV	40.59728	−115.379	Peromyscus maniculatus
UMNH	36546	United States, NV	40.59728	−115.379	Peromyscus maniculatus
UMNH	36547	United States, NV	40.62915	−115.367	Peromyscus maniculatus
UMNH	36548	United States, NV	40.59728	−115.379	Peromyscus maniculatus
UMNH	42480	United States, UT	39.9045	−113.37	Peromyscus maniculatus
UMNH	42481	United States, UT	39.9045	−113.37	Peromyscus maniculatus
UMNH	42482	United States, UT	39.9045	−113.37	Peromyscus maniculatus
UMNH	42483	United States, UT	39.88778	−113.413	Peromyscus maniculatus
UMNH	42484	United States, UT	39.88778	−113.413	Peromyscus maniculatus
UMNH	42485	United States, UT	39.88778	−113.413	Peromyscus maniculatus
UMNH	42486	United States, UT	39.88778	−113.413	Peromyscus maniculatus
UMNH	42487	United States, UT	39.88778	−113.413	Peromyscus maniculatus
UMNH	42488	United States, UT	39.88778	−113.413	Peromyscus maniculatus
UMNH	42489	United States, UT	39.88778	−113.413	Peromyscus maniculatus