Persistent chemicals from industrial processes, particularly perfluoroalkyl substances (PFAS), have become pervasive in the environment due to their persistence, long half-lives, and bioaccumulative properties. Used globally for their thermal resistance and repellence to water and oil PFAS have led to widespread environmental contamination. These compounds pose significant health risks with exposure through food, water, and dermal contact. Aquatic wildlife is particularly vulnerable as water bodies act as major transport and transformation mediums for PFAS. Their co-occurrence with microplastics may intensify the impact on aquatic species by influencing PFAS sorption and transport Despite progress in understanding the occurrence and fate of PFAS and microplastics in aquatic ecosystems, the toxicity of PFAS mixtures and their co-occurrence with other high-concern compounds remains poorly understood, especially over organisms’ life cycles.

Our study investigates the chronic toxicity of PFAS and microplastics on the sentinel species Daphnia, a species central to aquatic foodwebs and an ecotoxicology model. We examined the effects of perfluorooctane sulfonate (PFOS), perfluorooctanoic acid (PFOA), and polyethylene terephthalate microplastics (PET) both individually and in mixtures on Daphnia ecological endpoints. Unlike conventional studies, we used two Daphnia genotypes with distinct histories of chemical exposure. This approach revealed that PFAS and microplastics cause developmental failures, delayed sexual maturity and reduced somatic growth, with historical exposure to environmental pollution reducing tolerance to these persistent chemicals due to cumulative fitness costs. We also observed that the combined effect of the persistent chemicals analysed was 59% additive and 41% synergistic, whereas no antagonistic interactions were observed. The genotype-specific responses observed highlight the complex interplay between genetic background and pollutant exposure, emphasizing the importance of incorporating multiple genotypes in environmental risk assessments to more accurately predict the ecological impact of chemical pollutants.

Chronic toxicity of single chemicals and mixtures  

The toxicity of PFOS (70 ng/L), PFOA (7 ng/L), PET (50 mg/L) and their two-way and three-way mixtures were tested on the two Daphnia genotypes for the duration of their life cycles (until they released their second brood) (Fig. 1). To ensure a homogeneous PET solution, a stock solution of 5g/L was prepared and gently stirred on a magnetic stirrer before dosing the experimental vials at 50 mg/L.

Before the chronic toxicity test, we performed a qualitative assessment of Daphnia’s ingestion and egestion of microplastic over 72h. We exposed three clonal replicates of the two genotypes used subsequently in the chronic toxicity tests to 100mg/L of PET over 72h. Ingestion of microplastics was monitored every 24h in absence of feed. This experiment was run alongside a control in which Daphnia were fed algae. After 72h, Daphnia was transferred to clean medium every 24h, to visualise egestion of PET by microscope imaging. All images were captured using J software (Rueden et al., 2017). The medium was renewed every 24 hours and up to 72h after the last PET feed to estimate the gut retention time of PET and the time to full egestion.

Clonal replicates of the two genotypes were obtained from the stock collection of the University of Birmingham where they are maintained in standard laboratory conditions: 16:8-hr light–dark photoperiod; 0.8 mg L−1 Chlorella vulgaris fed weekly; ambient temperature: 10°C (Cuenca Cambronero et al., 2018). The growth medium used was borehole water, collected from a deep aquifer well and showing stable physico-chemical properties. Before commencing the exposures, clonal replicates of the two genotypes were acclimated in common garden conditions for three generations to the following conditions to reduce interference from maternal effects and synchronize reproduction: 16:8-hr light–dark photoperiod; 0.8 mg L−1 C. vulgaris fed daily; ambient temperature: 20± 2 °C. After three generations in these conditions, 24-hr-old randomly selected juveniles from the second or following broods were assigned to experimental conditions.

A total of 64 exposures were completed, including seven treatments and one control group for both genotypes and four clonal replicates per genotype and condition (Fig. 1). During the experiment, the growth medium (borehole water) was replenished every other day and spiked with the same concentration of chemical at each medium change to ensure a constant exposure throughout the experiment. Fitness-linked life history traits were measured during the experiment covering the Daphnia’s life cycle (until release of the second brood): age at maturity (first time the parthenogenetic eggs are released in the brood pouch); size at maturity (distance from the head to the base of the tail spine); fecundity (number of juveniles across the first two broods); interval between broods (days between the first and second brood) and mortality. Size was measured after the release of the second brood (end of the test) using ImageJ software (Rueden et al., 2017). 

Statistical analyses 

Univariate reaction norms were used to assess the effect of genotype (G), treatment (T -PFOS; PFOA; PET; PTE+PFOA; PTE+PFOS; PFOA+PFOS; PTE+PFOA+PFOS) and their interaction terms (G x T) on the five fitness-linked life history traits (age at maturity, size at maturity, fecundity, interval between broods and mortality) using a two-way ANOVA. A nested linear mixed effect model (LMMs) was used with clonal replicates as random effects nested within genotype using the “lmerTest” package in R version 4.3.1 (Kuznetsova et al., 2017). Before applying the LMMs model, the data were checked for normality using the Q-Q plots (residual vs. fitted value) (Zuur et al., 2010).   

The main effects of treatment, genotype and their interactions were visualised using ‘reaction norm’ plots. The mortality rate per treatment and genotype was estimated with the survival model fit using the psm function in the “rms” package in R version 3.6.0 (Harrell Jr, 2023). The day of mortality and mortality events were combined as a response variable while the term “genotype” was treated as a fixed effect. The mortality curves per generation were plotted with the survplot function from the rms package in R v.6.7.1 (Harrell Jr, 2023). A separate model was fitted to each treatment.

Multivariate effects on the overall fitness were calculated using multivariate statistics (MANOVA) by combining the life history traits (age at maturity, size at maturity, fecundity, and interval between broods) as response variable (y) and treatment and genotype as fixed terms (y ~ treatment * genotype). As per the ANOVA analysis of individual fitness-linked life history traits, a nested linear mixed effect model (LMMs) was used with clonal replicates as random effects nested within genotype using the “lmerTest” package in R version 4.3.1 (Kuznetsova et al., 2017). The multivariate reaction norms were visualized using phenotypic trajectory analysis (PTA) plots to describe the difference between control and treatments in terms of magnitude and direction of change following (Adams and Collyer, 2009).

 Mixture effects

In the common garden experiments, we investigated whether the chemical mixtures (two-way and three-way combinations) exhibited synergistic, antagonistic, or additive effects on Daphnia ecological endpoints (fitness-linked life history traits). To assess this, we used a null model of additivity (Cote et al., 2016), which contrasts the expected additive effect of two or more compounds (i.e., the sum of each compound individual effects) on ecological endpoints with the empirical observed effect of mixtures. These inferences were possible because the same set of genotypes were exposed to the single chemicals and mixtures.

The null addictive prediction of the joint effect of stressors combinations was calculated as: Emix=EA+EB+EC

where EA, EB and EC are the effects of single compounds and Emix is the mixtures effect. For each trait ‘y’ (age at maturity, size at maturity, fecundity, and interval between broods) standardized effect sizes were used corresponding to (y_treatment - y_control)/σ, where σ is the shared standard deviation of the pooled control and treatment trait values per treatment following (Cote et al., 2016). This analysis could not be done on mortality because individuals killed by one chemical cannot be killed by another simultaneously and, hence, the additive model cannot apply to this trait. The effect of mixtures was considered additive if the predicted joint effect ‘Emix’ was within the 95% confidence intervals of the observed effect from the multiple stressors and not significantly different from the null model. The effect of multiple stressors was antagonistic if the predicted joint effect was significantly smaller of the null model and synergistic if the predicted joint effect was significantly larger than the null model. We used Welch Two sample t-test to assess significant departure from the null model. This test compares the means of the predicted null model and the empirical observed values for each fitness-linked life-history trait.