Study Population
The Qanuilirpitaa? 2017 Nunavik Health Survey (Q2017) survey was conducted from August 19 to October 5, 2017, and covered all 14 Nunavik villages situated across three ecological regions. Accessing these isolated coastal villages was made possible through the Canadian Coast Guard Icebreaker, Amundsen, which facilitated the collection of data as participants were invited aboard the ship. The survey aimed to include all Nunavik permanent residents aged 16 years and above, and a structured stratified proportional model was employed for selecting respondents. A total of 1326 residents were successfully recruited for the study. Data collection encompassed comprehensive questionnaires, clinical measurements, as well as the collection of biological samples such as urine and blood.
The Q2017 survey adhered to the principles of OCAP® (Ownership, Control, Access, and Possession), ensuring respectful collaboration with various Nunavik organizations throughout the entire process. Inuit colleagues and partners and the Q2017 Nunavik Surveys Committee contributed to all discussions and result interpretation. The Q2017 Steering Committee, comprising members from diverse Nunavik organizations and the 14 municipalities, oversaw all research utilizing the dataset. Additionally, several co-authors played roles in co-designing, managing, and implementing the data collection process for Q2017. The Q2017 survey received ethical approval from the Comité d’éthique de la recherche du Centre Hospitalier Universitaire de Québec - Université Laval (#2016–2499 and #2016–2499-21).
Exposure biomarkers
The blood samples were collected by research nurses via venipuncture using K2-EDTA vacutainers. Subsequently, the collected samples were centrifuged at 2000×g for 10 minutes, and the resulting plasma was transferred into polypropylene tubes for storage at − 20°C until further analysis.
The analysis of perfluoroalkyl and polyfluoroalkyl substances (PFAS), including perfluorobutanoic acid (PFBA), perfluorohexanoic acid (PFHxA), perfluorooctanoic acid (PFOA), perfluorononanoic acid (PFNA), perfluorodecanoic acid (PFDA), perfluoroundecanoic acid (PFUnDA), perfluorobutane sulfonic acid (PFBS), perfluorohexane sulfonic acid (PFHxS), and perfluorooctane sulfonic acid (PFOS), was conducted at the Centre de Toxicologie du Québec (CTQ) of the Institut national de santé publique du Québec (INSPQ), which holds accreditation from the Canadian Association for Environmental Analytical Laboratories and ISO 17025.
For the PFAS analysis, plasma samples (100 µL) were enriched with labeled internal standards (PFBA-13C4, PFHxA-13C6, PFOA-13C4, PFNA-13C9, PFDA-13C9, PFUndA-13C7, PFHxS-13C3, and PFOS-13C4) and acidified using a 50% formic acid solution. Subsequently, solid phase extraction (SPE) was employed with SiliaPrepX WAX cartridges 100 mg/3 mL (Silicycle; Québec, Canada) to extract the analytes. The cartridges were initially washed with 5% NH4OH in methanol to remove contaminants and then conditioned with methanol and 2% formic acid before processing the samples. Following the loading of samples on the cartridges, a washing step with a 2% formic acid solution and methanol was performed, and the analytes were subsequently eluted using 3 mL of 5% NH4OH in methanol. The extracts were evaporated to dryness and then reconstituted in 1 mL of 5 mM ammonium acetate in 20% acetonitrile.
Analysis of the samples was carried out using Ultra Performance Liquid Chromatography (UPLC Waters Acquity) coupled with tandem mass spectrometry (MS/MS Waters Xevo TQ-S) (Waters, Milford, MA, USA) operating in the multiple reaction monitoring (MRM) mode with an electrospray ion source in the negative mode. The analytical column utilized was an ACE EXCEL C18-PFP 100 mm × 2.1 mm, 2.0 µm (ACE; Aberdeen, Scotland). The mobile phase consisted of an acetonitrile:H2O gradient (10:90) with 5 mM ammonium acetate, which was gradually changed to 100% acetonitrile with 5 mM ammonium acetate over a 7.0-minute run at a flow rate of 0.5 mL/min. This analytical method represents an improved version of the one previously described by Caron-Beaudoin et al. (2020), and both methods were cross-validated, demonstrating equivalent results. The limit of detection (LOD) for the perfluorinated compounds ranged between 0.01 and 0.07 µg/L. Intra-day precision varied between 2.0 and 3.8%, while inter-day precision varied between 3.8 and 5.0%, depending on the analytes. To monitor contamination, two laboratory blanks containing demineralized water were included in each analytic sequence. Polypropylene materials, such as tubes and microvials, were utilized to prevent the adsorption of PFAAs, and most of the other materials (glassware, pipette tips, etc.) were washed with methanol before use to minimize contamination. Additionally, an isolator column was inserted before the analytical column to control potential contaminants originating from the UPLC system.
To ensure the quality of the analyses, internal reference materials, such as the certified reference material SRM-1958 from the National Institute of Standards and Technology (NIST; Gaithersburg, MD), and in-house quality controls (QCs) for PFAS were employed. The overall quality and accuracy of the analytical method were monitored through participation in interlaboratory programs, including the AMAP External Quality Assessment Scheme (Centre de Toxicologie du Québec (CTQ), Institut National de Santé Publique du Québec (INSPQ), Québec, Canada) for the analytes PFHxA, PFOA, PFNA, PFDA, PFUnDA, PFHxS, and PFOS, as well as the German External Quality Assessment Scheme (G-EQUAS; Erlangen, Germany) for the analytes PFOA and PFOS.
Furthermore, the analysis of red blood cell fatty acid composition was conducted at the Laboratory of Nutritional Lipidomics of the University of Waterloo, Ontario, employing a Varian 3900 gas chromatograph equipped with a 15 m DB-FFAP capillary column (df = 0.10 µm) and a flame-ionisation detector. Concentrations of RBC eicosapentaenoic acid, docosahexaenoic acid, and docosapentaenoic acid were summed and expressed as a percentage of total fatty acids. This variable was further categorized into quartiles (n-3 PUFA). Serum vitamin D were determined using a MODULAR ANALYTICS e170 from Roche Diagnostics GmbH (Mannheim, Germany).
Whole blood total mercury (a surrogate for methylmercury exposure in fish and marine mammal eating populations) concentrations were measured using inductively coupled plasma mass spectrometry (ICP-MS) with a NexION® instrument from PerkinElmer (Waltham, MA, USA) and were analyzed at the Centre de toxicologie du Québec (Quebec, QC). Further details have been described elsewhere 30.
Outcome measurement
Respiratory health was measured via questionnaires (self-reported symptoms), spirometry measures, and medical records.
Self-reported symptoms included wheezing, chronic cough, chronic sputum, and breathlessness. Wheezing was defined as "wheezing or whistling in the chest at any time during the last 12 months", and is considered an indicator of asthma 31. There is no direct translation of the term "wheezing" in Inuktitut, as such this was translated to "laboured breathing" () or "whistling sound" (). The reliability of this indicator was previously tested against airway obstruction and reported an OR 2.26 (95% CI 1.46, 3.48) 11. The most commonly used epidemiological definition for chronic cough was adopted, defined as a "cough on most days for three months each year" 32. Chronic sputum was defined as “bringing up mucus on most days for three months each year, if the person usually brings up mucus from his chest when she does not have a cold”, and breathlessness was defined as “walking slower than people of the same age on the level because of breathlessness or having to stop for breath when walking at its own pace on the level”, corresponding to the second stage of the Modified Medical Research Council scale 33.
The main spirometry measures included in this analysis were the forced expiratory volume in the first second (FEV1) and forced vital capacity (FVC). The FEV1/FVC ratio below the lower limit of normal (LLN) (a continuous variable) was used to examine airway obstruction 34. No reference equation exists for the Inuit population, however the FEV1/FVC ratio shows minimal variation across ethnic groups 35. The LLN was calculated based on age and sex using the NHANES III reference equation for Caucasians 11,36. A second airway obstruction variable was calculated based on the fixed ratio recommended by GOLD (< 0.7) (a binary variable). It serves as an epidemiological definition for Chronic Obstructive Pulmonary Disease (COPD) in adults over 35 years old and can also be indicative of asthma, suboptimal lung development, or early lung decline 37.
Medical files (outpatient visits, emergency visits and hospitalizations) were reviewed by trained research nurses for physician diagnoses of asthma, COPD (including emphysema and chronic bronchitis), tuberculosis (latent or active/pulmonary) and hospitalization for respiratory infection during childhood (before the age of 5 years). Diagnosis of asthma or COPD in medical file are not necessarily supported by spirometry because it is not widely available in Nunavik.
Covariates
Participant sociodemographic characteristics were gathered through questionnaires. The covariates of interest included age (years), education level (< 9th Grade level, at least some high school, at least some college), personal income (<$20,000, $20,000-$59,999, >$60,000), marital status (married or common law, single), smoking status (never, former, current smoker < 15 pack-years, current smoker ≥ 15 pack-years), second-hand smoking exposure (< 1/month, 1/week – 1/month, nearly everyday), marijuana use (rarely, 1–3/month, ≥ 1/week), food security (secure, moderately insecure, severely insecure), house crowding (≤ 1 person-per-room, > 1 person-per-room), vitamin D concentrations, n-3 PUFA in red blood cell quartiles (a marker of seafood and marine mammal consumption), fruit and vegetable consumption (< 5 times/day, ≥5times/day) 38 , waist circumference quartiles, and blood total mercury concentrations (Supplementary Fig. 1).
Dietary factors associated with respiratory outcomes and PFAS were included in the models. Among the main exposure sources of vitamin D and n-3 PUFA in red blood cells are marine mammals and seafood 27,28,39,40 which are also the main exposure sources of PFAS.
Further details on these covariates (calculations and questionnaire responses) have been described elsewhere 11,29. In place of body mass index (BMI), we used waist circumference as a proxy for body size. This decision was influenced by the fact that BMI tends to overestimate the prevalence of obesity in Inuit populations 41. Some recent evidence suggests a link between exposure to mercury and respiratory health 42–44. Given that elevated mercury exposure is common in Nunavik due to the high mercury concentrations found in some marine mammals and fish consumed by Inuit 27,30, we included whole blood total mercury levels as a potential confounding factor.
Statistical analysis
Multiple imputation using fully conditional specification was utilized to impute missing data with plausible values. The technique involved multivariate linear regression for imputing continuous variables, logistic regression for imputing binary variables and multinomial logistic regression for imputing categorical variables 45. Predictions for each variable were made with the inclusion of all other variables used in the study, including socioeconomic, lifestyle and respiratory predictors. Given that 57% of participants had missing values for at least one variable (typically a few variables), a total of 60 complete datasets were imputed. Each imputed dataset underwent separate analysis, and the results were subsequently pooled using the "mianalyze" SAS procedure. Descriptive data from non-imputed and imputed datasets showed high similarity. Notably, respiratory outcomes were not imputed and participants with missing outcome data (n = 38 for wheezing, n = 39 for chronic cough, n = 208 for airway obstruction) were excluded from the analysis.
Plasma levels of all PFAS congeners were log2-transformed. Geometric means and their corresponding 95% confidence intervals (CIs) were calculated for exposure variables, and proportions were calculated for the categorical covariate and outcome variables. Covariate and outcome variables were compared by sex using t-tests. The percentage of samples with detectable PFBA, PFHxA, and PFBS concentrations were low with > 80% of values below the LOD and were not included in further analyses. Study demographics were explored in the non-imputed dataset to compare the distributions in imputed and non-imputed datasets. We also calculated the distributions in an imputed and weighted dataset. Survey weights were included to account for the complex stratified sampling strategy and the non-response rate 9,46.
Binomial generalized linear models was conducted to explore associations between individual PFAS congeners and respiratory outcomes wheezing, airway obstruction, and asthma. The study focused on wheezing, chronic cough, airway obstruction, and asthma because these were considered to be respiratory symptoms of concern in the Nunavik context 11. However, we also examined associations between PFAS and chronic sputum, breathlessness, and having any respiratory symptom in supplementary analyses. Beta coefficients were exponentiated to calculate the odds ratio (and associated 95% CI) for a doubling of the exposure. Additionally, multivariate linear regression models were used to examine the associations between PFAS congeners and FEV1, FVC and the FEV1/FVC ratio. Beta coefficients were exponentiated to calculate the unit outcome change for every doubling of the exposure.
Covariates were included in two stages. The first adjusted model controlled for sociodemographic and some lifestyle variables: sex, age, personal income, marital status, smoking status, second-hand smoking, marijuana use, waist circumference, overcrowding, and food security. The second adjusted models were further adjusted for n-3 PUFA red blood cell quartiles, vitamin D concentrations, fruit/vegetable intake, and blood mercury concentrations to account for other nutrient and environmental contaminants detected in PFAS exposure sources. No issues of multicollinearity emerged after checking correlation matrices and model variance inflation factors (VIF) values.
We also examined the modifying effect of nutrition on the associations between PFAS and respiratory outcomes, namely n-3 PUFA in red blood cell quartiles and vitamin D deficiency (defined as vitamin D levels < 30 ng/mL) (described above). An interaction term was included in all models between individual PFAS congeners and each of the nutritional variables. If a significant interaction was observed (p-value < 0.05), we stratified the models by the nutritional categories.
Outcomes associated with PFAS from single chemical models were further explored using mixture analyses. The association between PFAS congeners and respiratory outcomes was modeled by Bayesian Kernel Machine Regression (BKMR), which allows for correlation, non-linearity and interactions between exposures. The BKMR model relies on a function of PFAS congeners to estimate the exposure-response surface used to compute contrasts of interest and inferences. For example, the overall effect of the exposure mixture can be obtained by comparing the point estimates of the exposure-response surface when all mixture components are set to a particular quantile versus all set at a reference level (50th quantile). Furthermore, exposure-response functions can be obtained to visualize non-linear relationships and to investigate interactions between mixture components by estimating the bivariate exposure-response function by varying level (25th, 50th and 75th quantiles) of another component. The bkmr R package 47 was used to estimate the covariate-adjusted exposure-response surface and to calculate credible intervals (CI) for inferences. Missing data were imputed in 60 datasets and a separate BKMR model with 10,000 iterations in the Markov chain Monte Carlo algorithm was fit for each dataset. The first 5,000 iterations were discarded and 20% of remaining iterations were used to compute posterior means and variances. Estimated of each dataset were combined using Rubin’s method 48. The theoretical background of the BKMR method and implementation with this package has been described in detail elsewhere 47,49.
Sensitivity analyses
According to the directed acyclic diagram (DAG) created (Supplementary Fig. 1), adjusting for participants with recent respiratory infections or surgery may lead to bias by collider. We did, however, conduct a sensitivity analysis excluding participants with emphysema or recent/past history of tuberculosis. Furthermore, asthma is sometimes misdiagnosed with chronic obstructive pulmonary disease (COPD) in older adults. As such, we conducted another analysis expanding the asthma definition to include COPD diagnosis, self-reported asthma, self-reported asthmatic bronchitis, or self-reported allergic bronchitis. We also replaced the current smoking status variable with cotinine levels and found no differences in the results.
The cross-sectional nature of the study limits our ability to establish temporality, particularly with regards to asthma. While asthma diagnosis in adulthood is common 50, most people with asthma are diagnosed as children 51. To better ascertain temporality, we conducted additional sensitivity analyses by restricting the analysis of PFAS and asthma to those aged 16–20 years. Given the persistent nature of PFAS, we assumed that PFAS levels remained relatively stable throughout childhood. This assumption is strengthened in the Nunavik context since the main exposure source of PFAS is consumption of country (traditional) foods and is reflected by ties to culture, so we assumed that the consumption of country versus market foods remained relatively stable in childhood as well.
All regression model analyses were performed using SAS 9.4 (SAS Institute, Cary, NC, USA) and BKMR models were performed using R 4.3.1 software.