Bone Density and Genomic Analysis Unfold Cold Adaptation Mechanisms of Extinct Inhabitants of Tierra Del Fuego

The Fuegians, extinct inhabitants of Tierra del Fuego, were an exemplary case of a cold-adapted population capable of living in extreme weather conditions without any adequate clothing, however the mechanisms of their extraordinary resistance to cold remain enigmatic. Brown adipose tissue (BAT) plays a crucial role in this kind of adaptation besides having a protective role on the detrimental effect of low temperatures on bone structure. Skeletal remains of 12 adult Fuegians were analyzed for bone mineral density and structure. We show that, despite the unfavorable climate, bone mineral density of Fuegians was close to that seen in BAT-expressing modern humans of temperate zones. Furthermore, we report signicant differences between Fuegians and other cold-adapted populations in the frequency of the Homeobox protein Hox-C4 (HOXC4) rs190771160 variant, a gene involved in BAT differentiation, whose identied variant is predicted to upregulate HOXC4 expression. Greater BAT accumulation might therefore explain the Fuegians extreme cold-resistance and the protection of their skeletons against major cold-related damage. Our results increase our understanding of how ecological challenges have been important drivers of genetic and environmental factors interactions.


Introduction
Subsequently to the dispersal of Homo sapiens out of Africa during the Paleolithic [1], human populations required adaptations to diverse climatic conditions [2] achieved through morphological and cultural adjustments and metabolic mutations [3], such as brown adipose tissue (BAT) upregulation in cold climates [4]. The interest in understanding the mechanisms of cold adaptation comes both from obvious anthropological and evolutionary implications, and from theories that attribute to populations with higher BAT activity an increased resistance to obesity and diabetes [5]. When Europeans reached Tierra del Fuego in 1520, it was inhabited by different populations (mostly named Yamana, Alacaluf and Ona) generally grouped under the term Fuegians, today extinct [6]. What is known mostly comes from navigators and ethnographers' reports (e.g., Bove G. 1884; Gusinde 1966) [7,8]. All were nomadic hunter-gatherers, and they had an exceptional cold resistance: despite the extremely harsh climate, they did not use any closed clothing, and open res and animal fat smeared on bare skin appeared inadequate means of protection to those describing them [6].
Cultural and physiological mechanisms implicated in cold adaptation include clothing and shelter, peripheral vasoconstriction, muscle and subcutaneous fat-mediated insulation, lower surface area-to-mass ratio [9], thyroid hormone levels variation [10,11], high-fat and -protein diet [10,12], BAT over-activity, and white adipose tissue browning [9,13,14]. In humans, BAT increases upon cold exposure [15,16], and it generates heat via upregulation of uncoupling protein-1 (UCP1), increasing basal metabolic rate (BMR) [17]. Interestingly, in 1960 Hammel studied nine Fuegians shortly before the extinction of this population, showing that their BMR was about 160% of that reported for the cold exposed "average white man" (p. 24) [18]. This report raised the suspicion that the Fuegians might have had increased BAT activity at the root of their amazing cold adaptation.
Beyond the interest in the role of BAT in thermoregulation and as a target for the prevention and treatment of obesity and type 2 diabetes, coldactivated BAT is also positively related to bone mineral density (BMD) [19] and femoral cross-sectional area (CSA) [20]. BAT volume parallels BMD [21] and predicts femoral CSA and cortical bone area (CBA) [22]. Moreover, BAT-impaired cold-exposed mice show more bone loss than WT via β2 adrenergic receptor activation [23], and UCP1 knock-out (KO) mice under permanent cold stress have lower bone formation and mineralization compared to WT [24]. These altered bone phenotypes are not seen in thermoneutral UCP-1 KO and WT mice, indicating that BAT may prevent an otherwise massive cold-induced bone loss [24].
To test the hypothesis that BAT could have contributed to their exceptional cold adaptation, we evaluated bone morphology and density of the skeletal remains of the collection of Fuegians preserved in the Museum of Anthropology 'G. Sergi' of the Sapienza University of Rome (Italy), and we analyzed genomes of this extinct population available in public repositories for variants in genes involved in BAT modulation.

Results And Discussion
To con rm BAT and bone physiology interconnection, widely reported in animal models but less established in human beings, we rst explored the correlation between BAT and femoral cortical thickness in a living population composed of 34 BAT expressing individuals from a previously described large cohort [15]. In line with previous reports, BAT volume and activity were directly related to femoral CBA (Table S1). Moreover, BAT volume independently predicted CBA after height and muscle area (Table S2). After con rming this association in the living population, we analyzed Fuegians' bones for comparison.
The Fuegians skeleton series preserved in Rome is composed of 14 complete skeletons; it represents one of the largest collections outside Argentina and is in an exceptionally good state of preservation. Twelve adults were included, and their bones were analyzed (Table S3). Fuegians femoral CBA and CSA were signi cantly lower than those of BAT positive and negative living humans, whereas the endocortical bone area (EndA) was comparable (Fig. 1A). A dual-energy X-ray-absorptiometry (DXA) of the Fuegians femur and lumbar spine remains showed no different BMD compared with a subpopulation of 12 matched BAT expressing subjects living in temperate zones and selected among the cohort described above (Fig. 1B).
The observation that, despite their exposure to environmental conditions unusually detrimental for bone, the Fuegians showed similar BMD compared to BAT expressing modern subjects living in temperate zones supports the hypothesis that a high BAT expression might have protected their skeletons from even more negative consequences. Yet, the con icting nding of a lower CBA and CSA led us to further test whether the Fuegians had an abundance of BAT through a genomic approach.
We therefore analyzed available genomes of this population, and we obtained reliable results for ve of the 11 ancient Fuegian genomes retrieved [25]. Coding and non-coding variants in a panel of 28 genes involved in BAT functional pathways [13] were identi ed in the Fuegian genomes (DatasetS1). All variants found in the Fuegians were searched in a control population of 14 individuals including 10 Siberians, 2 Athabascan and 2 Greenlanders [25]. These individuals, all living in cold areas of the globe (modern Alaska, Greenland and Siberia), were selected for the less striking cold adaptation leading to different lifestyle and habits compared to Fuegians, such as the use of warm clothing.
Interestingly, a coding variant previously associated with cold adaptation [13], the cAMP-dependent protein kinase type II-beta regulatory subunit (PRKAR2B) rs75385144, was identi ed in four out of ve Fuegians (Table 1). This candidate variant has a very low frequency in the African population (0.07) compared with European, American and Asian populations (0.30-0.40) [13], however it was found in a frequency comparable to Fuegians in our control population (11 out of 14 individuals, p = 0.95). These results suggest that this variant might be implicated in cold adaptation, but not in different ways across these populations.
Moreover, a statistically signi cant difference in frequency between Fuegians and the control population emerged for two non-coding variants, PR domain containing 16 (PRDM16) rs2493270 and Homeobox protein Hox-C4 (HOXC4) rs190771160, both genes being involved in BAT differentiation [13,26]  PRDM16 is present in brown adipocytes and is able to upregulate UCP-1 expression, thus increasing the activity of BAT. Its overexpression is also associated with browning in mouse models [27]. Its variant rs2493270 was found in three out of ve Fuegians and was not identi ed in controls (p = 0.002) ( Table 1, Fig. 2). In gnomAD, the variant is reported with a very low frequency in the Finnish population, and in slightly higher frequencies in Latinos, East Asians and Africans. The variant is an A to G substitution within the 3' UTR region of the PRDM16 gene, predicted as likely benign by in silico tools, although a role in modifying genomic interactions and miRNA binding sites cannot be excluded with possible functional consequences on beige and brown adipocytes.
HOXC4 rs190771160 was found in four out of ve Fuegians and only one out of 14 controls (p = 0.0015) ( Table 1, Fig. 2). This variant is reported with a very low frequency in all modern populations whose genomes are deposited in gnomAD database. The variant is a C to T substitution within a highly phylogenetic conserved intronic/promoter region of HOXC4. Intriguingly, all of the interrogated in silico tools consistently predicted a very likely causative effect in affecting histone modi cations, mainly tri-methylation at lysine 27 of histone H3 (H3K27me3).
Moreover, the transcription factor USF1 binding site seems to be affected by the variant as predicted by the same in silico tools (Table 1). Enhanced expression of HOXC4 has been linked to the accumulation of fat droplets within brown adipocytes [26]. Anticipated to affect H3K27me3 and USF1 binding, both inhibiting transcription [28,29], the presence of rs190771160 is likely to upregulate HOXC4 expression, leading to enhanced BAT proliferation in those carrying the variant.
Overall, both the genomic data and bone analysis of the Fuegians suggest that this population might have been characterized by large BAT depots. However, the possibility that the Fuegians could have gained their resistance to cold through other functional pathways cannot be de nitely ruled out. Likewise, we cannot exclude that factors such as muscle mass, whose positive correlation with bone mass is established, may have contributed to their better-than-expected BMDs. However, even if this was the case, our hypothesis could be maintained since the myokines derived from skeletal muscle after contraction or cold exposure (i.e. irisin) seem to have profound effects on enhancing bone mass and reinforcing the BAT phenotype [30][31][32][33].
In conclusion, although based on a small sample, our ndings may con rm the hypothesis that the extraordinary cold-adaptation of the Fuegians, as noted by the nineteenth-century explorers, was due to an exceptional BAT accumulation.

Materials And Methods
Study populations

Fuegian remains
The skeletal remains of 13 adult Yàmanas from Tierra del Fuego were selected from the collection of the Museum of Anthropology "G. Sergi" of the Sapienza University of Rome (Italy), which includes fteen complete or almost complete skeletons (13 adults, 1 juvenile, 1 infant) in good state of preservation. Great part of the collection was recovered in 1883-1884 by the Italian explorer Giacomo Bove during one of his expeditions in South America, and then acquired by the University of Rome (1886) through a donation for the rising museum curated by Giuseppe Sergi [34].
The selection criteria included the availability of femur and lumbar vertebras and estimated age of 25 years old or older as proof of achievement of adult bone structure.

Living population
We retrospectively selected a population of 217 patients expressing F-FDG BAT among a total of 6454 patients who underwent 8004 consecutive 18F-FDG PET/CT scans from January 2007 to June 2010 at the Istituto Nazionale dei Tumori Regina Elena (Rome, Italy, 41.81°N 12.45°W). The modality of BAT expression detection is described elsewhere [15]. Among the 217 patients with BAT depots, individuals younger than 25 years old or who underwent PET/CT scans for malignant diseases and received the last treatment for malignant diseases within one year before the scan were excluded, as where subjects whose diagnosis and/or treatment timings were unknown. We also excluded those with diabetes, renal failure, and under steroids or osteoporosis treatment at the time of PET/CT. Moreover, a gender, age and BMI matched control population of BAT negative subjects was selected among the same population following the same exclusion criteria stated above. To control for BMI difference, given the reported lean body mass of the Fuegians, we only selected patients with a BMI ≤25 kg/m 2 .
Among the patients enrolled in the present study, a subpopulation was subsequently selected to be age, BMI and gender matched to those estimated for the fuegian skeletal remains adopting the following criteria: available Dual energy X-ray absorptiometry (DXA) scan performed within 1 year of the PET scan with a Hologic scanner.
Positron emission tomography/computed tomography (PET/CT) of living humans PET/CT scans were performed using a Biograph 16 High Rez PET/CT scanner (Siemens AG, Munich, Germany). Patients came from the temperate metropolitan area of Rome and were instructed to fast overnight for at least 12 h before the scan and to abstain from carbohydrates and very fatty foods consumption, nicotine, caffeine, or alcohol intake for the preceding 24 hours. Room temperature water was allowed at all times. Since their arrival to the hospital, the patients were in an air-conditioned environment at about 22°C. After the injection of 18 Fuorodeoxyglucose (FDG) the patients rested at 24°C for 1 hour before undergoing the PET/CT scan. Data on gender, BMI, age, plasma glucose level and malignancy status (active: PET/CT scan positive for malignancy, not active= PET/CT scan negative for malignancy) were obtained for all patients. Parameters measured from PET/CT scans included, in addition to BAT characteristics: cross-sectional area (CSA, mm2), cortical bone area (CBA, mm2), endocortical bone area (EndA, mm2), muscle area (MA, mm2) taken at proximal-shaft diaphysis, 1 cm distal to the lesser trochanter. These were calculated with the use of ImageJ 1.52i, Wayne Rasband, National Institutes of Health, USA.

Computed tomography (CT) of Fuegians femurs
For each Fuegian skeleton, the right femur was analyzed using a commercial CT scanner (Revolution, GE Healthcare). Before the scan, femurs were oriented in standardized anatomical planes. Diaphyseal structural properties were calculated from DICOM les with the use of ImageJ 1.52i (Wayne Rasband, National Institutes of Health, USA). Parameters measured included: cross-sectional area (CSA, mm2), cortical bone area (CBA, mm2), and endocortical bone area (EndA, mm2) taken at proximal-shaft diaphysis, 1 cm distal to the lesser trochanter.
DXA of femur and lumbar spine of Fuegians Lumbar spine (LS, L1-L4), femoral neck (FN) and total hip (TH) BMD of Fuegians femur and lumbar vertebrae remains were measured through DXA (QDR Discovery Acclaim, Hologic Inc., Waltham, MA, USA). Before DXA scanning, femurs of Fuegians were oriented in standardized anatomical planes and the lumbar vertebrae placed in a rack that allowed them to be aligned, using rice as a soft tissue proxy [35].

Genome analysis
Alignment data from Raghavan et al., 2015 [25] were downloaded from http://www.cbs.dtu.dk/suppl/NativeAmerican/. Data corresponding to 8 Fuegians and 14 distinct Native Americans individuals exposed to cold temperatures were identi ed and extracted. Read statistics per sample and alignment methods are shown in the Supplementary Materials for Raghavan et al. [25].
Since the individual Fuegians libraries have a variable endogenous content (from 0.6% to 23.8% of the total number of reads mapped to the human genome) and a low sequencing depth (average depth from 0.003 to 1.7) [25], samples with less than 5 million of mapper reads were removed prior to further analyses.  [13], were selected and further investigated (S1).
Candidate variants were further characterized based on possible functional effect as predicted by in silico analysis using MutationTaster [37]and RegulationSpotter [38].These are online applications performing several in silico tests on both DNA and protein level ultimately estimating the impact of the identi ed variant, such as functional consequences of synonymous or intronic mutations up to amino acid substitutions, deletion and insertion of sequences, or variants within the intron-exon border.

Ethics statement
Data relative to the study cohort of living humans was collected from a database previously used to report BAT prevalence in central Italy [15]. The project had been approved by the Institutional Review Board of the Dipartimento di Prevenzione e Diagnostica Oncologica, UOC di Medicina Nucleare, Istituto Nazionale Tumori Regina Elena, Roma. The data were analyzed anonymously. Written informed consent was given by the patients or from the next of kin or caretakers for their information to be stored in the hospital database and used for research. All methods were carried out in accordance with relevant guidelines and regulations.

Statistics
Statistical analysis was performed using SPSS 25.0 (SPSS, Inc., Chicago, IL). Data are expressed as mean ± standard deviation for normally distributed variables. Variables were tested for normality of distribution using Shapiro-Wilk test. Variables that were not normally distributed were log-transformed. Comparisons between groups were performed using the paired Student's t-test or Wilcoxon test as appropriate. Relationships between variables were measured by Pearson's correlation coe cient. Stepwise regression modeling was used to determine predictors of CBA. A chi-square test was performed to compare variant frequencies between Fuegians, Siberians/North Americans, and other modern populations from gnomAD database (https://gnomad.broadinstitute.org/). An α error of 0.05 was considered the threshold for statistical signi cance Declarations