Chemicals

The following PFASs were tested in the present study: perfluorooctanoic acid (PFOA), perfluorononanoic acid (PFNA), perfluorohexane sulfonate (PFHxS), and perfluorooctane sulfonate (PFOS). In addition, the effects of the FOXO1 inhibitor AS1842856 were studied. All stocks were prepared in 100% dimethyl sulfoxide (DMSO HybriMax, Sigma-Aldrich). More information about suppliers, purity, catalog numbers, CAS numbers and maximum concentrations tested in the present study is presented in Supplementary Table 1.

Namalwa cell culture

The human Burkitt’s lymphoma cell line Namalwa was obtained from Sigma-Aldrich (Zwijndrecht, The Netherlands) and maintained in RPMI1640 (Gibco, Thermo Fisher Scientific, Waltham, MA) supplemented with 10% heat-inactivated fetal bovine serum (Gibco), 1% Sodium Pyruvate (Sigma–Aldrich), 1% NEAA (Gibco) and 1% penicillin–streptomycin (Sigma–Aldrich). Once the Namalwa cells were growing at a stable growth rate (approximately after 10 days), the amount of heat-inactivated fetal bovine serum was reduced to 2% to limit binding of PFASs to proteins to allow significant cell exposure. Namalwa cells were subcultured twice a week, each time diluted to 0.5 × 106 viable cells/mL.

Cell viability studies

The effects of the 4 PFASs and AS1842856 on Namalwa cell viability were determined using the WST-1 assay. This assay determines the conversion of the tetrazolium salt WST-1 (4-[3-(4-iodophenyl)-2-(4-nitrophenyl)-2H-5-tetrazolio]-1,3-benzene disulfonate) to formazan by metabolically active cells. Namalwa cells were cultured in 96-well plates (1 × 10^6 cells/ml) and exposed to increasing concentrations up to 100 µM for PFOA, PFOS, PFNA and PFHxS for 48 h or to increasing concentration up to 0.1 µM AS1842856 for 48 h. After exposure, WST-1 solution (Sigma-Aldrich) was added to the cell culture medium (1:10 dilution). After 1 h incubation in an incubator (humidified atmosphere with 5% CO2 at 37 °C), the plate was shaken at 1000 rpm for 1 min, and absorbance at 450 nm was measured (background absorbance at 630 nm was subtracted) using a Synergy HT Microplate Reader (BioTek, Winooski, VT). Cell viability upon PFAS and AS1842856 treatments was expressed as percentage of the cell viability of the solvent control. DMSO concentration was 0.1% in solvent control and in all treatment conditions.

Namalwa exposure for gene expression analysis

For gene expression studies, Namalwa cells were seeded in 24-well plates (Corning, NY; 0.5–1 × 10^6 cells per well in 500 µL). Test chemicals were first diluted from a 1000-fold concentrated stock to a twofold concentrated stock solution in assay medium and subsequently twofold diluted upon the addition to the Namalwa cell suspension, providing a final DMSO concentration of 0.1%. In each experiment a solvent control (0.1% DMSO) was included. For the RNA seq study, Namalwa cells were exposed to 100 µM PFOA for 48 h. For RT-qPCR studies, Namalwa cells were exposed to either the highest non-cytotoxic concentration of 100 µM for PFOA, PFOS and PFHxS, and 33 µM for PFNA for 6, 24 and 48 h or to a concentration range up to 100 µM for PFOA, PFOS and PFHxS, and up to 33 µM for PFNA for 48 h. In a follow-up study, Namalwa cells were exposed to 0.0001, 0.001 and 0.01 µM AS1842856 for 6, 24 and 48 h. After exposure, effects of the PFASs and AS1842856 on expression of selected genes were assessed.

RNA library preparations and RNA seq

To obtain insight into the PFOA-induced gene expression changes, Namalwa cells were exposed to 100 µM PFOA or solvent control for in total 48 h. After exposure, medium including cells were collected and centrifuged at 200 × g for 7 min to obtain a cell pellet. The Namalwa cell pellets were subsequently lysed in cell lysis buffer (RLT) and total RNA was isolated and purified using the RNeasy Mini Kit (Qiagen). Total RNA was quantified using Qubit (Life Technologies) and RNA integrity was analyzed using Agilent 2100 Bioanalyzer total RNA Pico chip (Agilent Technologies). Subsequent RNA Library preparations and RNA sequencing was performed at Genomics Facility of Wageningen University and Research, Business Unit Bioscience. Approximately, 1 µg total RNA was used for RNA library preparation using TruSeq Stranded mRNA Sample Prep kit (Illumina). In short, after polyA based mRNA selection, RNA were further processed including subsequent fragmentation, first and second strand cDNA synthesis, adapter ligation and final library amplification resulting in RNA seq libraries including unique dual indexes, all following manufacturer’s protocol. Final libraries were eluted in 30 µl elution buffer followed by library quality assessment using a Fragment Analyzer (Agilent Technologies) and quantified by Qubit fluorescence measurements (Invitrogen, Life Technologies).

Prepared libraries were pooled in an equimolar manner and combined with other indexed libraries for sequencing on an Illumina NovaSeq 6000 system. Final sequencing was done using an S2 and S4 type flow cell, both with XP loading workflow and settings specific for 2 × 150 nt paired end reads plus dual indexes reads. All steps for sequencing were carried out according to manufacturer’s protocol. Demultiplexing of reads per sample by corresponding indexes was done using bcl2fastq v2.20.0.422 (Illumina Inc, San Diego, CA, USA).

Processing of RNA sequencing reads

The RNA seq reads were used to quantify transcript abundances. To this end, the tool Cutadapt (version 1.16) (Martin 2011) was used to trim adapters from the reads and HISAT2 (version 2.2.1) (Kim et al. 2019) was used to map the reads to the GRCm38.13 human genome assembly-based transcriptome sequences as annotated by the GENCODE consortium (release M40). HISAT2 output was converted and sorted by chromosomal position using Samtools (version 1.9) (Danecek et al. 2021). RSeQC (v3.0.1) (Wang et al. 2012) and PRINSEQ (v0.20.4) (Schmieder and Edwards 2011) were used for quality control. HTSeq (version 0.11.2) (Anders et al. 2015) was used to count reads in transcripts using gene-level quantification. Differential gene expression was determined using the package limma (version 3.50.3) (Ritchie et al. 2015) utilizing the obtained scaled gene-level counts. Briefly, before statistical analyses, nonspecific filtering of the count table was performed to increase detection power, based on the requirement that a gene should have an expression level greater than around 10 counts in at least three samples. Differences in library size were adjusted by the trimmed mean of M-values normalization method, implemented in the package edgeR (version 3.36.0) (Robinson et al. 2009; McCarthy et al. 2012; Chen et al. 2016). Counts were transformed to log2 (cpm) values and associated precision weights, and entered into the limma analysis pipeline. Differentially expressed genes were identified using generalized linear models that incorporate empirical Bayes methods.

Ingenuity pathway analysis

Gene lists containing gene identifiers (Ensembl Gene ID), and corresponding log2 fold changes and p-values were uploaded to Ingenuity Pathway Analysis (IPA) software (Qiagen, Redwood City, CA, USA). Input criteria were a log2 fold change of above 0.5 and a p value below 0.01. To interpret biological meaning of differentially expressed genes, gene sets were analyzed using the Canonical Pathways module. Statistically overrepresented pathways were identified by Fisher’s exact test (p < 0.01).

RT-qPCR

After exposure of Namalwa cells to PFASs, medium including cells were collected and centrifuged at 200 × g for 7 min to obtain a cell pellet. The Namalwa cell pellets were lysed in cell lysis buffer (RLT) and total RNA was extracted using the RNeasy Mini Kit (Qiagen, Venlo, The Netherlands). Subsequently, 500 ng RNA was used to synthesize cDNA using the iScript cDNA synthesis kit (Bio-Rad Laboratories, Veenendaal, The Netherlands). Changes in gene expression were determined by RT-qPCR on a CFX384 real-time PCR detection system (Bio-Rad Laboratories) using SensiMix (Bioline; GC Biotech, Alphen aan den Rijn, The Netherlands). The PCR conditions consisted of an initial denaturation of 95 °C for 10 min, followed by 40 cycles of denaturation at 95 °C for 10 s and annealing extension at 60 °C for 15 s. Primer sequences were taken from the Harvard PrimerBank and ordered from Eurogentec (Liège, Belgium). Sequences of the used primers are listed in Table 1. Relative gene expression was quantified with the standard curve method, using a standard curve generated from a serial dilution of pooled sample cDNA, and subsequently normalized to glyceraldehyde 3-phosphate dehydrogenase (GAPDH) gene expression. Gene expression upon PFAS treatments was expressed as fold change compared to the gene expression measured for the solvent control. The concentration–response data were subjected to BMD analysis using PROAST software as described as follows.

Table 1 Primer sequences used for RT-qPCRBenchmark dose (BMD) analysis of RT-qPCR data using PROAST

RT-qPCR data were used for concentration–response modeling with BMD analysis software PROAST version 70.5 (National Institute for Public Health and the Environment 2018) in R (version 4.2.0). PROAST is particularly applied for modeling of in vivo (dose–response) data, providing information on the BMD. In the present work, PROAST software was used for the analysis of in vitro (concentration–response) data, thereby providing information on the benchmark concentration (BMC). Data of all PFASs were analyzed simultaneously to ensure the parallel curves required to derive RPFs. Tab-delimited text files containing data on concentration, experiment number, mean effect, standard deviation, and sample size (number of replicates) were made and analyzed as continuous (summary) data. For PROAST analysis, average gene expression values of triplicates with standard deviation of two independent studies were used (Supplementary file 2). Independent studies were assigned as covariates. Then, the exponential model:

with y denoting the response and x the concentration was applied. The parameters a, b, c, and d describe the response at dose 0 (background value), the potency of the PFAS, maximum fold change in response compared with background response (upper or lower plateau), and steepness of the curve (on a log-dose scale), respectively. BMC values were determined for a benchmark response of 10% (BMR10) for the RT-qPCR data. The model with the lowest Akaike information criterion (AIC) was chosen to determine RPF values including 90% confidence intervals. PFOA was used as the index chemical. Of note, PROAST definitions are CES (critical effect size) and CED (critical effect dose), which are the same as BMR and BMC, respectively.

Statistical analysis

Data are presented as mean ± SD. A one-way ANOVA followed by Dunnett’s post hoc multiple comparison test was used for comparisons between Namalwa cells exposed to various concentrations of AS1842856 and solvent control. p < 0.05 was considered as statistically significant. Prism software (version 9.2.0; Graphpad, San Diego, CA) was used for statistical analysis.