No evidence of inbreeding depression in a Tasmanian devil insurance population despite significant variation in inbreeding

Inbreeding depression occurs when inbred individuals experience reduced fitness as a result of reduced genome-wide heterozygosity. The Tasmanian devil faces extinction due to a contagious cancer, devil facial tumour disease (DFTD). An insurance metapopulation was established in 2006 to ensure the survival of the species and to be used as a source population for re-wilding and genetic rescue. The emergence of DFTD and the rapid decline of wild devil populations have rendered the species at risk of inbreeding depression. We used 33 microsatellite loci to (1) reconstruct a pedigree for the insurance population and (2) estimate genome-wide heterozygosity for 200 individuals. Using heterozygosity-fitness correlations, we investigated the effect of heterozygosity on six diverse fitness measures (ulna length, asymmetry, weight-at-weaning, testes volume, reproductive success and survival). Despite statistically significant evidence of variation in individual inbreeding in this population, we found no associations between inbreeding and any of our six fitness measurements. We propose that the benign environment in captivity may decrease the intensity of inbreeding depression, relative to the stressful conditions in the wild. Future work will need to measure fitness of released animals to facilitate translation of this data to the broader conservation management of the species in its native range.

and not inbred, as information from the wild is generally unavailable. Rudnick & Lacy (2008) explored the impact of founder relationship assumptions on the success of mean kinship management using simulations. The long-term consequences (over 5 generations) of assuming founders to be unrelated were small when compared to a model of known founder relationships and known inbreeding values. However, short-term consequences (3 to 5 generations) of an initial increase in inbreeding values were exposed when founders were assumed to be unrelated compared to a model of known founder relationships and inbreeding 7 . This finding highlights the importance of understanding the recent ancestry of founders when commencing an insurance population, particularly those which will be used for rewilding and genetic rescues soon after instigation, such as in the case of the Tasmanian devil insurance population.
Genetic diversity of the Tasmanian devil is low 8,9 . Tasmanian devils were previously found on mainland Australia approximately 5,000 years ago, but are now restricted to the island state of Tasmania 10 . At least three population crashes, resulting in genetic bottlenecks, have been documented in Tasmania since devils were isolated on the island approximately 10,000 years ago 11 . Most recently, this species has suffered an extreme population decline (~85%) 12 due to the emergence of two transmissible cancer lines, Devil Facial Tumour 1 (DFT1) and Devil Facial Tumour 2 (DFT2), first documented in 1996 and 2015 respectively 13,14 . These severe population declines have rendered the Tasmanian devil at risk of further genetic diversity loss and inbreeding depression.
Inbreeding depression is the decrease in individual fitness that occurs as a result of increased homozygosity via inbreeding. The scale of the Tasmanian devil insurance population, and access to genetic samples from the entire population, provides a unique opportunity to explore the effects of inbreeding on fitness in an intensively managed, newly established population. Inbreeding depression is well documented in both captive and wild populations across numerous taxa (e.g. refs 15 and 16) and can vary widely both among and within species, for multiple reasons (e.g. founding population size and diversity 17,18 or available resources 16 ). This can create a dilemma for captive breeding programs, as the general strategy for management is to minimise mean kinship. However as each species program and life-history varies an understanding of the relationship between inbreeding, life-history and management practice may be needed to minimise the negative effects of inbreeding.
The detection and evaluation of inbreeding, and inbreeding depression, can be achieved through two different approaches: pedigree-derived inbreeding coefficients, or molecular-based inbreeding estimates. The relationship between these two metrics and the relationship they independently have with individual-fitness measures have differing strengths and weakenesses [19][20][21][22] . Wright (1992) provided a method to calculate a coefficient of inbreeding, F, using pedigree information. This approach relies on two assumptions; (1) founding individuals of the pedigree are unrelated, and (2) the pedigree information is an accurate representation of genetic relationships and degree of shared ancestry. If portions of the pedigree are unknown, or the relationship between recent founders is unknown, pedigree-derived inbreeding values inaccurately represent true variation in genetic relationships. Molecular-based inbreeding values offer an alternative to pedigree-based values, by providing a multi-locus heterozygosity score intended to be representative of genome-wide inbreeding levels 22 . This technique operates predominately under the expectation that increased inbreeding will increase the probability that an individual is homozygous for deleterious recessive alleles, leading to inbreeding depression 20,21 . This approach however can be limited by sampling variance, such that multi-locus heterozygosity at the markers used may not be representative of genome-wide heterozygosity. Nevertheless, molecular measures of inbreeding enable heterozygosity-fitness correlations (HFCs) to be used as a tool for detecting inbreeding depression in individuals with limited or unknown pedigree information.
HFCs may be based on a wide range of direct and indirect measures of fitness, generally categorised as direct measures of fitness, including life-history traits (survival and reproductive success), versus indirect measures of fitness, including physiological (e.g. parasitic load) and morphological (e.g. body weight) traits 19 . On average, slightly stronger HFCs are detected when using direct fitness measures 19 . Significant HFCs are still detected using indirect fitness correlates [23][24][25] , and these measures are often more easily attainable than direct measures, such as survival of long-lived species, or reproductive success in wild populations. In the presence of inbreeding depression, HFCs may not be detected unanimously across all measures of fitness (e.g. refs 26 and 27). Thus, studying both direct and indirect mechanisms by which inbreeding depression affects various traits in a single species is necessary. Identifying these traits can help conservation research more broadly.
Ultimately, the detection of HFCs can only take place when there is variance in inbreeding among the study population 20,28 . Variance in inbreeding within a population can be detected as identity disequilibrium. Using a multi-locus approach, covariance in heterozygosity across markers can be used to measure identity disequilibrium and quantify the variance in inbreeding 21, 29 . In endangered species, particularly those in captive breeding programs derived from restricted founding populations, the likelihood of identity disequilibrium is high.
In this study, we investigated whether a captive population (part of the Tasmanian devil insurance metapopulation) showed signs of inbreeding depression across six measures of fitness using both pedigree-based and molecular-based approaches. We developed microsatellite markers to enable pedigree reconstruction of devil juveniles produced in group-breeding enclosures (multiple males and females housed together), and examined whether the pedigree or molecular data are informative of inbreeding depression at six fitness traits, including direct and indirect measures of fitness, in captivity.

Methods
Study Species and Data Collection. This study was conducted at a Tasmanian devil breeding facility located in NSW, Australia. The facility consists of 12 managed environmental enclosures (MEEs) and eight intensive enclosures built over 25 hectares. In MEEs, Tasmanian devils are housed at a density of approximately 2.67 devils/hectare, consisting of eight 3-hectare enclosures holding 4:4 (male:female) devils and four 2.5-hectare enclosures holding 3:4 devils 30 . Intensive enclosures generally house 1-2 individuals, either for welfare reasons or breeding purposes. Breeding recommendations are issued annually and not every devil is provided with a breeding opportunity every year. The breeding facility was established in January 2011 with 30 devils (11 wild founders 3 , 2 hand-raised orphans, and 17 captive born devils from the existing zoo insurance population 2 ). A limited number of transfers into the facility occurred in subsequent breeding seasons (fewer than 10 devils per year). However, 26 of the transferred individuals were related (litter mates), comprising 10 litters.
Tasmanian devil females give birth to multiple offspring (20)(21)(22)(23)(24)(25)(26)(27)(28)(29)(30), but have only four teats in their pouch and as a result a maximum of four offspring can survive per female per breeding season 31 . Mixed paternity litters are possible. The breeding season occurs between March to June 32 . Sampling for this study occurred during the annual catch up event over a three-week period during November and December in 2014 (under approval from University of Sydney Ethics Committee #2014/550). Following approved ethical standards, Tasmanian devils were restrained in a hessian bag and ulna length, asymmetry and testes volume were recorded (full measurement methods are below). As part of regular husbandry practices (following approved ethical standards), reproductive outputs and weight at weaning were recorded for all Tasmanian devils from 2011 to 2014.
Due to the nature of fitness measurements collected, six datasets were identified for statistical analysis. Table 1 describes the fitness measurement collected for each dataset and the corresponding sample size. In total, 200 Tasmanian devils were included in this study.
Phenotypic Traits. Six phenotypic traits were measured (Table 1) Each fitness trait was analysed independently for inbreeding depression in this study. Along with the above-mentioned fitness measurements, each individual had the following data recorded: age, sex, year, pen total (the number of individuals housed within the same enclosure) and pen ID. This data was obtained from studbook information 38 and annual housing reports from the breeding institution.

Marker Development, DNA Extractions and Molecular Measurements. DNA samples had been
previously extracted from ear biopsies and blood samples collected as part of the management of the Tasmanian devil insurance population 39 with further samples collected in 2014 under the direction of the ZAA, on behalf of the STDP, for management purposes. Samples were extracted using a modified phenol-chloroform extraction technique 40 and stored at −20 °C. All individuals were genotyped for 33 microsatellite markers (Table 2).
At the commencement of this study there were only 11 putatively-neutral microsatellite markers characterised in devil 41 (Table 2). In order to assess HFCs, additional microsatellite marker development was required. Devil-specific microsatellite makers were developed using the Tasmanian devil genome accessible on Ensembl 42 . Microsatellite repeats were identified throughout the genome using RepeatMasker 43 . Genomic data processed in RepeatMasker was a minimum of 10 kb away from coding DNA and sourced evenly across all six Tasmanian devil autosomes. All candidate microsatellites were a minimum of 10 kb apart, to minimise the possibility of linkage disequilibrium. In total, 44 candidate microsatellites were chosen for polymorphism screening. Primer pairs were designed and optimised using Oligo 7 (Molecular Biology Insights). Loci were initially screened for variation using 12 Tasmanian devils, randomly selected from the insurance population. For the initial screening process, forward primers were labelled with a universal 6-FAM CAG Tag (CAGTCGGGCGTCATCA) (following 44 ). PCRs were performed following the standard protocol for the Qiagen Type-it Microsatellite kit, with modified primer concentrations of 0.06 µM for tagged primers and 0.6 µM for untagged primers and the CAG tag. Thermocycling conditions followed a protocol of 95 °C for 5 min, followed by 30 cycles of 95 °C for 30 s, 60 °C for 90 s and 72 °C for 30 s with a final extension of 60 °C for 30 min. Twenty-two polymorphic microsatellites were identified ( Table 2): five on autosome 1, six on autosome 2, seven on autosome 4, one on autosome 5 and three on autosome 6 (5 and 6  For genotyping all animals included in this study, the forward primer of each marker was fluorescently labelled ( Table 2). PCRs were performed following the standard protocol for the Qiagen Type-it Microsatellite kit, with the modification that the reaction size totalled 10 µL. Fragments were separated on an ABi3300 using GeneScan ™ 500 LIZ ™ as size standard. Nine multiplexes were developed using Multiplex Manager 45 ( Internal relatedness (IR) was calculated as a measure of individual multilocus heterozygosity, using the Rhh 51 package in R 52 . The IR metric gives higher weighting to rare alleles in the population, making this method ideal for studying inbreeding 53 . Homozygous individuals for a rare allele will have a higher IR value than homozygous individuals for common alleles, due to the increased chance that rare allele homozygosity is the product of an inbred mating 53 . There are numerous methods to calculate multilocus heterozygosity, a commonly used alternative being standardised multilocus heterozygosity. In our dataset, the correlation between IR and standardised heterozygosity was 0.908. As the correlation was high and IR is most suitable for studying inbreeding 19 , we selected IR for subsequent analyses.
Pedigree Construction. Parentage (paternity and maternity) of offspring produced in group pens was assigned using CERVUS 54 , using a combination of molecular data and annual enclosure records of the location of each Tasmanian devil. Paternal and maternal candidates included all sexually mature Tasmanian devils present at the institution during the breeding season for each offspring cohort. A total of 118 joeys were born during the years 2011-2014 (2011 n = 24, 2012 n = 38, 2013 n = 29 and 2014 n = 27); all of which were genotyped in this study. All inferred parent-pair combinations were checked against annual reports and only accepted if individuals had been present in the same enclosure during the reproductive season. Critical LOD scores were calculated for each cohort by simulation of 10,000 offspring at an error rate of 1% (proportion of loci mistyped = 0.01). LOD confidence calculations were set at strict confidence = 95% and relaxed confidence = 80%. All candidate parents were analysed using joint LOD scores to assign parentage. All genetically confirmed parents (i.e. known parents) had 0 mismatches with their offspring. All remaining potential parent pairs had >3 mismatches, and were therefore confidently excluded as potential parents for a given offspring. Known parentage data (genetically assigned parents) was entered into the studbook and pedigree analysis undertaken using the PMx software 55 .
Inbreeding Analysis. Identity disequilibrium can arise after, among other reasons, inbred matings, and result in heterozygosity and homozygosity correlations throughout the genome, within a population (reviewed in ref. 21). To assess whether a chosen panel of markers (number and type) inform individual inbreeding status the markers should be representative of variation in identity-by-descent (IBD) 29,21 . Pedigrees may fail to detect variance in inbreeding due to unknown founder relationships, while molecular measures of IBD can expose inbreeding where it previously would not have been detected. Statistical evidence of identity disequilibrium provides support that a chosen marker panel measures variation in identity-by-descent and therefore may be useful for interpreting heterozygosity-fitness correlations under inbreeding 21 . Two analyses are commonly used to quantify the degree of identity disequilibrium, heterozygosity-heterozygosity correlations (HHCs) 56 and the g 2 statistic 57 . Both methods involve quantifying the degree of covariance in heterozygosity among markers to test for identity disequilibrium; we used both approaches here. We calculated g 2 with the RMES software 57 , using 1000 permutations. The R package Rhh was used to calculate HHC using 1000 Monte Carlo simulations.
Statistical Analysis. Our general approach used generalised linear mixed models (GLMMs) to assess the impact of internal relatedness (our predictor of interest) on five of our measures of fitness: ulna length, asymmetry,  reproductive success, weight at weaning and testes volume (survival was modelled differently, see below). Mixed models were fitted using the R package lme4 58 , where our indirect fitness traits of ulna length, asymmetry, weight at weaning and testes volume were fitted with lmer (Gaussian response variables) and our direct fitness traits (reproductive success) were fitted with glmer (binomial response variables). Fixed predictors in each global model included IR, as well as a unique set of fixed and random factors as relevant to the biology of each trait and each dataset (see below; full model specifications provided at Supplementary Table S2). Each global model was standardised to facilitate comparison of parameter estimates across models using the arm package 59 in R. Model selection proceeded under information theory (following 60 ). Briefly, we created a complete subset of models using the MuMIn package 61 in R, submodels were ranked using the Akaike's information criterion (corrected for sample size; AIC C ), and model averaging was used to take all models that fell within 2 AIC C of the highest ranked model. Inference was based on standardised effect sizes, their standard errors, and the relative importance scores for each included parameter 60 . This method can result in the exclusion of IR as a predictor from the final model, if its explanatory power is low. Under inbreeding depression, we predicted IR to be negatively correlated with testes volume, asymmetry, ulna length, weight at weaning, reproductive success and survival.

Indirect Fitness Traits.
In addition to the above general specifications, our indirect fitness trait models also included pen ID as a categorical random factor in each model (to account for between-pen variation). Our weight at weaning model also included "maternal IR" and "paternal IR" (in addition to the individual offspring IR) as continuous fixed predictors in the global model, due to the strong potential influence of parental effects on early juvenile traits. Specifically, the maternal investment provided to pouch young, and the molecular effects of both maternal and paternal internal relatedness on joey development. In the latter model, individuals were only included if all three IR values were available (joey IR, maternal IR and paternal IR). Global models for morphological traits also included sex (43:48 M:F), age and body weight (mean = 7.31 kg, SD = 1.62, n = 91) where appropriate (Supplementary Table S2). The age range of both male and females in this dataset spanned juveniles (1 year of age) to 5-year-old individuals. Tasmanian devils do not breed past the age of 5 32 and pre-DFTD in the wild survived to 5 to 6 years of age 62 . They are known to live up to 9 years of age in captivity 63 .

Direct Fitness Traits. Reproductive success for both females and males was modelled as a binomial process.
For females, the number of events (successes) was the number of offspring produced by a female, and the number of trials was the total number of offspring biologically possible per breeding attempt (four potential joeys). For males, the number of events (successes) was the number of offspring sired by a male, and the number of trials was the total number of offspring sired by all individuals in the pen. Representing reproductive success as a proportional value in this way standardises across enclosures of differing population densities, sex ratios, and overall productivity (especially for males). Enclosure reproductive success was overdispersed for both sexes, whereby typically one male sired the majority of offspring in an enclosure, and approximately ~50% of females successfully reproduced in an enclosure. To account for overdispersion in our models we therefore added a residual parameter for each individual. In addition to IR, the male and female reproductive success models also included age as a continuous variable and pen ID and year as categorical random factors. For both male and female devils the captive environment does not truly offer an opportunity for free mate choice, because males have limited opportunity to maximise their reproductive success, and females may choose mates from only a select group of males, which can vary in heterozygosity across enclosures and across years. We investigated whether female mate choice was for the most heterozygous male per enclosure (each male is compared to the remaining males in their enclosure), or alternatively, if female mate choice had a heterozygosity threshold (regardless of heterozygosity ranking per enclosure, all individuals with a certain level of homozygosity are generally not chosen as mates). To investigate this question, we ran the same model previously described for male reproductive success, however we standardised male internal relatedness within enclosure (referred to herein as "s.IR"). This was calculated as the difference between a males individual IR and the mean IR for all males within the enclosure (i.e. s.IR = IR − IR enclosure mean ) (following 64 ).

Survival Analysis.
A Cox proportional hazard model was used to investigate the effect of inbreeding on survival, with the R package survival 65 . This analysis allows for the inclusion of animals that outlive the study period. Age of death was recorded for all Tasmanian devils that died during the study period. Tasmanian devils were censored if they were still alive at 1 March 2015 (study end date), with their age recorded. The covariates of IR and sex were fitted as predictors of survival. The null expectation of hazard rate (HR) is equal to 1, where the hazard of dying increases with estimates >1, and survival probability increases with estimates <1.

Results
We found no evidence of null alleles and no significant deviations from Hardy-Weinberg (Table 2). Missing data was low: out of the 200 Tasmanian devils included in this study, ≥95% were genotyped for ≥90% of the 33 markers. Scoring errors based on regenotyping of approximately 16.5% of individuals were minimal (0.2%). The mean number of alleles per locus was 4.00 (range 2-9, n = 33 loci, 200 devils). Mean observed heterozygosity was 0.447 ± (SD) 0.167, and mean expected heterozygosity was 0.459 ± (SD) 0.170 (Table 2).
Scientific RepoRts | 7: 1830 | DOI:10.1038/s41598-017-02000-y Pedigree reconstruction. All 118 joeys born during 2011-2014 were successfully assigned parentage using joint LOD scores for all candidate parents. Cohabitation was found to be true of the most likely parent-pair combination for 100% of the joeys examined. Incorporating full pedigree information from the insurance population, we observed that 98.6% of individuals in our study population had a pedigree derived inbreeding coefficient of 0.00, with only three greater than 0 (0.03, 0.11 and 0.25) (mean = 0.002 ± SD 0.020, n = 200). As variance is required in the degree of inbreeding to investigate associations among pedigree-based inbreeding and fitness, the pedigree data was uninformative and not analysed further 20 .
Heterozygosity-fitness correlations. After model averaging, internal relatedness (IR) appeared in the top model or a model within 2 AIC C of the top model for testes volume and male reproductive success (Table 3). Using standardised IR for males gave broadly similar results to absolute IR (Table 3). Standard errors were large for IR in models where it was included (confidence intervals encompassed zero) and relative importance values of IR were low (0.27-0.71; Table 3). For the remaining fitness components (weight at weaning, ulna length, asymmetry and female reproductive success) IR was not included in the top model, nor a model within 2 AIC C of the top model. Of the 97 individuals with weaning weight data only 79 had IR values for joey IR, maternal IR and paternal IR, the remaining were excluded from the model. Missing data comprised of two joey IR values, nine maternal IR values and seven paternal IR values. As missing data was found across all three IR categories, we proceeded to keep our data restricted to joeys that had IR values across all three categories.
In our final models, fitness traits were largely influenced by the sex and age of individuals (these parameters had a relative importance value of 1 for several responses) ( Table 3). We observed effects of sex on ulna length (males were larger, while controlling for age) and weight at weaning (males were larger, while controlling for age)    (Table 3). Unsurprisingly, older individuals showed longer ulna measurements (Table 3). Female reproductive fitness declined with age, but male reproductive fitness did not ( Table 3). The survival hazard rate analysis indicated that there is an imprecise increased risk of the event "death" occurring, with increasing IR. This result was not statistically significant (the confidence interval for this analysis includes the null hypothesis, 1 Table 4).

Discussion
We successfully assigned all offspring (using molecular data) to their respective parents and reconstructed a pedigree for a component of the Tasmanian devil insurance population. We investigated the effect of inbreeding on six fitness measures in this population and surprisingly, found no evidence for inbreeding depression. While the pedigree based values showed negligible evidence of inbreeding (98.6% of inbreeding coefficients = 0.00, with only three individuals having an F > 0.00), our HHC and g 2 analyses based on microsatellite genotyping revealed statistically significant variation in inbreeding in the Tasmanian devil insurance population. Inbreeding coefficients derived from our molecular data were consistent with trapping proximity (see map 3 ) which supports suggestions by Hogg et al. 3 that some original founders in the insurance population may have been closely related individuals. Understanding these potential relationships is beyond the scope of this study but is being addressed by others within our research group. In comparison to a simulation-based study 66 that evaluated the effectiveness of g 2 for estimating the variation in inbreeding, the significant g 2 value we report here fell within the range indicated as sufficient to detect correlations between survival and inbreeding 66 . Our results highlight the importance of using both pedigree and molecular tools for population management. A molecular-based approach can provide valuable information in regards to inbreeding that may not necessarily be reflected by a pedigree-based approach (e.g. ref. 67) or when a pedigree is unavailable. Molecular tools can also provide valuable information about individual relationships, such as parentage, which can be used to strengthen pedigree based management and breeding strategies going forward (e.g. ref. 68).
We used a range of life history and morphological traits in an attempt to maximise the possibility of detecting inbreeding depression in the Tasmanian devil, if it is occurring. Indirect fitness measures (both morphometric and physiological) generally experience stabilizing selection pressure toward an optimum, which can create a non-linear relationship between the fitness measurement and heterozygosity 19,21 . In contrast, direct fitness measures are generally under directional selection and involve a greater number of loci than morphometric and physiological measures 21, 29 . These two contrasting forms of selection can influence the detectability and strength of HFCs, where there is a greater likelihood of detecting HFCs in direct fitness measurements.
None of the direct or indirect fitness measures we recorded showed a strong association with IR in captive Tasmanian devils. There are several explanations for this pattern: (1) this population of Tasmanian devils exhibits little variation in inbreeding (e.g. all individuals have high, or low, levels of inbreeding with little-to-no variance), (2) variation in inbreeding is high, but is not adequately captured by the markers used in this study, or (3) inbreeding has only a weak effect on the fitness traits we measured, in this population. We have shown evidence that there is variance in inbreeding, as our g 2 and HHC results were both significant, indicating our genetic dataset does reflect inter-individual variation in inbreeding and thus has the potential to detect inbreeding depression if it occurs (negating arguments 1 and 2, above). Our results may reflect the potential for environmental effects in captivity to lessen inbreeding depression, i.e. reduce the variation in fitness between more/less inbred individuals, and our dataset may have insufficient power to detect such a weak effect.
Given that our molecular dataset has the capacity to detect inbreeding depression if it occurs, we feel that environmental effects may be a plausible explanation for our results in the Tasmanian devil. The effects of inbreeding depression can be decreased in benign environments relative to stressful environments 69,70 . For example, environmental variation, which may act as a stressor, can differ between a species' wild environment (generally more variable/stressful) and its captive environment (generally more consistent/benign). Armbruster & Reed (2005) conducted a meta-analysis of 34 studies investigating inbreeding-environment interactions and found an average fitness reduction of 69% in populations exposed to stressful environments relative to benign environments. Thus, inbreeding depression may be weaker in captive conditions, such as the population included in the current study, and may explain the minimal effect we observed between internal relatedness and trait variation 71 . Management practices for the institution in this study involves quarterly health checks (dental examination, wound treatment and parasitic treatment), a feeding regime that is directly proportional to the number of individuals housed in each enclosure, and an assortment of differing enclosure sizes with varying ratios of males and females (breeding facility manager, pers. comm.). Longevity, resource competition and reproductive competition can all differ in captivity relative to the wild.
A particular environmental mechanism that may impact our results is resource excess, which may mask any signs of inbreeding depression that would be detected under stressful conditions 72 . In general, parental investment can be influenced by resource availability, indirectly impacting reproductive success (e.g. refs 37 and 73). Tasmanian devils give birth to highly altricial young that undergo an intense period of growth in the mother's pouch 74 . From February to September (breeding season and pouch young), Tasmanian devils at the breeding facility in this study are fed excess (25% increase) resources to accommodate for pouch young investment (breeding facility manager, pers. comm.). Minimising resource competition among females in this captive population could mask the effects of competitive resource investments, explaining our observation that other measures, such as age, are better predictors of reproductive success than IR. This suggests that there is room for between-individual variation when the effect is large enough (the effect of biological aging) as opposed to small effect predictors (such as internal relatedness). Contrasting to females, resource investment plays a lesser role in male reproductive success, as male marsupials have negligible direct input to their offspring after fertilisation 74 . Tasmanian devil courtship behaviour is not well documented, but females do display both aggressive and submissive behaviours towards potential mates 75 . Forced copulation due to the larger size of male Tasmanian devils (average weight 8.46 kg [current study]) relative to female Tasmanian devils (average weight 7.27 kg [current study]) may result in males acquiring reproductive partners for reasons other than their own IR (e.g. body weight). Furthermore, the majority of Tasmanian devils in the current study were of optimal weight due to resource excess. As such, both low IR males and high IR males have equal opportunity for obtaining recourses, resulting in equal opportunity of acquiring a reproductive partner. For the Tasmanian devil, the full impact of inbreeding depression may only be fully realised under wild conditions. Future Directions and Management Implications. We reconstructed a pedigree and found minimal variation in inbreeding coefficients, whereas our molecular data indicated significant variation in underlying identity by descent (inbreeding). Our results suggest that some of the founders of the insurance population may have been closely related, as pairings informed by the shallow pedigree have apparently resulted in some mating among relatives. Our molecular observations of inbreeding, in the absence of pedigree inbreeding, highlight the need for future conservation programs to ensure that the relationships among founding individuals from small, fragmented populations are well understood prior to the commencement of a large-scale insurance population. Incorporating such data will ensure that any potential inbreeding is minimised and accounted for in future breeding management decisions. These considerations are of particular consequence to insurance populations that engage in re-wilding and genetic rescue programs early after establishment, such as the Tasmanian devil insurance population.
Until recently, obtaining an accurate genome-wide estimate of identity-by-descent (IBD) with molecular technology was difficult. However up-to-date research reports that molecular measures of inbreeding can be as accurate, and in some cases, more accurate than pedigree estimates 22,76 . Between-sibling variation in IBD is not reflected in pedigree-derived inbreeding values (which represent only 'expected' values), but can be detected using a sufficiently large panel of molecular markers. The integration of pedigree and molecular assessments will benefit captive population management (with respect to revealing founder relationships), and can be further applied to wild populations in the absence of any pedigree data.
Despite significant variation in inbreeding, we found no evidence for inbreeding depression in the Tasmanian devil; we propose that the benign environment in captivity may attenuate inbreeding depression. Tasmanian devils from the insurance population are already being released to both wild and semi-wild areas (peninsula and island sites) 2 and we plan to measure their fitness following release. Future work will enable us to test the hypothesis that highly inbred Tasmanian devils experience a fitness disadvantage when released into the wild.