Ethics
This is a multicenter, secondary analysis of a prospectively recruited longitudinal cohort study enrolling consecutive patients with suspected sepsis. All patients were enrolled under local ethical board approval. Informed written consent was obtained upon enrollment from the patient or their legal representative. The Clinical Research Ethics Board (REB) of the University of British Columbia (UBC) provided ethics approval for all sequencing and bioinformatics studies, carried out in a manner blinded to patient identity (approval number REB#H20-02441, REB#H17-01208). Patients recruited and enrolled at Unity Health Toronto were included in accordance with protocol approved by the St. Michael’s Hospital ethics board (REB#: 20-078). Patients’ data were extracted from the in-hospital electronic medical records, de-identified, and assigned random identification numbers which were used throughout the project. All experiments performed at the NRC involving human samples were approved by the NRC’s Ethics Board (NRC REB 2021-57) and experiments were performed according to NRC’s policies governing human subjects that follow applicable research guidelines compliant with the laws in the province of Québec.
Sample collection, RNA isolation and cDNA conversion
Patient samples were collected in Pax Gene tubes and total RNA was isolated using standard protocol for Qiagen RNAeasy mini kit (# 74104). RNA was assessed first using NanoDrop One spectrophotometer (Thermo Scientific) and A260/A280 values were between 1.8 and 2.2, with typical yields in the range of 6–8 µg total. RNA Integrity Number (RIN) was determined using the Agilent 2100 Bioanalyzer (Agilent Technologies). Following the standard Nanochip protocol, samples with RIN values > 7.0 were used for conversion to cDNA. Input volumes for reverse transcription were calculated using the concentration from the bioanalyzer (~500 ng total was used per sample) and a High Capacity cDNA Reverse Transcription Kit (Applied Biosystems # 4368814) was used following standard protocols. For RNA-Seq, results have been published33 including blood collection, RNA extraction and downstream processing. Accession numbers are included in Supplementary Table 1.
Discovery dataset
The whole transcriptome (RNA-seq) data from 586 whole blood samples from different countries and continents comprised our patient cohort. 514 samples were collected and used for discovery analyses (i.e., the discovery cohort). The remaining 72 samples were secondary samples collected from 72 individuals, which were excluded from discovery analyses (to prevent same-individual artifacts) and used as a validation cohort. The sepsis severity associated with the discovery cohort (514) was based on the SOFA score of the patient at 24 h after the first sample collection: 271 samples with SOFA ≥ 2 were sepsis, and 243 samples with SOFA < 2 were non-sepsis.
Sepsis signature and housekeeping gene selection
We tested 99 cellular reprogramming (CR) genes as potential sepsis markers34. We used DESeq2 to perform differential gene expression analyses and chose the genes that had the highest up-regulation (positive fold-change) in high severity samples. We also estimated the predictive accuracy of each CR gene by setting the sensitivity to 75%. We picked six genes (RETN, S100A8, MCEMP1, S100A12, CYP1B1 and HK3) that had the best results in both analyses for the Sepset model.
We analyzed 2833 housekeeping genes38 in our discovery cohort of 514 samples to set a baseline for RNA quantity and sequencing depth. We selected housekeeping genes (HKG) with high and consistent expression across all samples based on mean and variance. We then examined the expression variance of the top 20 HKG candidates across key clinical factors such as age group, gender, sepsis severity, patient location, mortality, etc. The two housekeeping genes (PTP4A2 and CHTOP) with the lowest variances were used to set a baseline for the SepsetER model.
ML algorithm construction and testing
Our own published RNA-Seq data from 873 patient samples33, was used for feature (gene) reduction using ML. An additional 1241 transcriptomes from patients were used for testing the derived signature. Three major groups of datasets were used for biomarker development – training, validation and testing. The discovery data set (N = 586) was first tested by 10X cross validation and randomly divided into a training (90% of the samples) used for the construction of the models (10,000+ models) and a test dataset (10% of the samples) to assess the best model (Supplementary Fig. 2). We trained 18 different machine learning algorithms on the transcriptomic profiles of the 514 discovery cohort samples. The algorithms were: K-Nearest Neighbors (KNN), Ridge Regression (RR), Lasso regression (LR), Elastic Net (EN), Partial Least Square (PLS), Linear Discriminant Analysis (LDA), Regularized Discriminant Analysis (RDA), Quadratic Discriminant Analysis (QDA), Bayesian Generalized Linear Model (BL), Naïve Bayes (NB), Support Vector Machines (SVM), Decision Tree (DT), Random Forest (RF), Adaptive Boosting (AB), Stochastic Gradient Boosting Model (GBM), Extreme Gradient Boosting (XGB), Neural Network (NN), and Multilayer Perceptron (MLP). We tested each model with different parameters and chose the best one based on the AUC-ROC using 10-fold cross validation, repeated 10 times. We then validated the performance of each model with the additional 72 validation samples (that were not in the training dataset).
We tested the Sepset model with multiple methods. We used various sepsis transcriptome datasets (from microarray and RNA-seq platforms) with over 3000 sepsis and healthy samples to evaluate the SepsetER sepsis classification model. We also trained other published sepsis gene-signatures with our training dataset and compared them with Sepset. The Sepset model, using the Extreme Gradient Boosting (XGB) algorithm, performed better than all other signatures, with a median AUC-ROC of 0.85 in all testing datasets.
Design of primers and probes
The expression of 6 top genes of interest was assessed based on the selection of the highest fold changes with respect to severity of disease in the ICU cohort. These genes (and their amplicon sizes) are HK3 (108 bp), RETN (78 bp), S100A12 (122 bp), S100A8 (122 bp), MCEMP1 (131 bp), CYP1B1 (114 bp). The housekeeping genes were also selected based on stable expression: PTP4A2 (138 bp) and CHTOP (113 bp).
Primers (IDT) for these genes have been designed to span the exon-exon junction to avoid amplification of genomic DNA and to cover the different isoforms. The amplicons’ sizes range between 78 bp and 138 bp, as stated above. The probes were synthesized (IDT) with either FAM, HEX, ROX or Cy5 fluorescent labels and a ZEN/3’ Iowa Black FQ (IABkFQ) double quencher, when possible, to reduce background noise.
Two fourplex reactions were designed to include 3 genes of interest and 1 HKG for normalization. As such one reaction targeted: CYP1B1, MCEMP1, S100A12 and PTP4A2, and the other: HK3, RETN, S100A8 and CHTOP. The sequences of the primers and probes are provided and described in the supplementary Table 5.
Specificity and sensitivity of the primers and probes were first assessed by performing qPCR standard curves of the individual targeted genes from human universal cDNA (P/N 637223, Clontech/TaKaRa Bio) and comparing the efficiency with the multiplex reaction. The results obtained were then used to design the multiplex ddPCR reactions in order to ensure appropriate amplification of differentially expressed genes and avoid amplification bias.
Commercial duplex ddPCR
For optimal results, recommendations made in the Droplet Digital PCR Applications Guide (Bio-Rad Bulletin 6407) were followed. We used equal concentrations of cDNA for droplet generation following the protocol for ddPCR supermix (Cat # 1863026, Bio-Rad). Briefly, a 22 µL reaction set-up consists of 2X supermix, 20X probes (a duplex of FAM and HEX), equal concentrations of the patient samples (500 pg), and RNase-free water. The bulk solution (in a 96-well plate) is applied to the AutoDG (automated droplet generator) where the solution is partitioned into 10,000 individual water-in-oil droplets. The 96 well plate is foil sealed and put into the C1000 thermal cycler (Bio-Rad) where the individual droplets are subjected to the following conditions: 10 min at 95 °C, 40 cycles of 30 s at 94 °C and 1 min at 60 °C, followed by 10 min at 98 °C and a 4 °C hold. Subsequently, the droplets were read in the QX200 Droplet Reader using FAM and HEX channel readout in the QuantaSoft software. After data acquisition, the QC of the samples was assessed (ensuring equal droplet numbers generated) and samples were selected in the well selector tool under the Analyze tab. Samples were all manually thresholded using the values from probe alone readout and confirmed in 2D tracings of the duplexed reaction. Samples were tested in duplicate. The concentration reported is “copies/ng DNA” of the final 1X ddPCR reaction.
Microfluidic device fabrication
Sample preparation cartridge
Microfluidic channels and reservoirs were carved into a block (50 mm × 100 mm × 6 mm) of Zeonor 1060 R (Zeon Chemicals) using precision machining (Q350 CNC Mill; Menig Automation). The machined polymer piece was cleaned with isopropanol (Sigma-Aldrich) and dried with a stream of nitrogen gas. The microfluidic circuit was sealed using adhesive film (ARclear 93495, 40 μm in thickness; Adhesive Research) applied on a polycarbonate sheet (#85585K103; 250 μm in thickness; McMaster-Carr).
Detection cartridge
The microfluidic circuit was fabricated in polydimethylsiloxane (Sylgard 184; Dow Corning) using replica molding. A multi-level SU-8/silicon master mold was made by sequential photolithographic patterning of multiple layers (10, 30 and 50 µm in thickness) of SU-8 photoresist (GM1060 and GM1070; Gersteltec) spin-coated onto a 6″ silicon wafer (Silicon Quest International) in conjunction with flood exposure at 365 nm (Hg i-line) through a chrome/quartz glass photomask (Photronics) using an EVG 6200 mask alignment system (EV Group). SU-8 resist was developed in propylene glycol monomethyl ether acetate (# 484431 Sigma-Aldrich) for several minutes, followed by rinsing with isopropanol (#I9030, Sigma-Aldrich) and drying with a stream of nitrogen gas. Bake steps were performed on a programmable hot plate (HS40A; Torrey Pines Scientific) using recommended time and temperature settings. The liquid pre-polymers of PDMS were mixed at a ratio of 10:1 (w/w) elastomer base/curing agent, poured onto the SU-8/silicon master mold, and cured at 85 °C for 1 h. The cured PDMS replica was bonded to a glass substrate following oxygen plasma activation (HI RF power, 900 mTorr for 30 s; Harrick Plasma).
Microfluidic assay implementation
Total RNA extraction from whole blood
Total RNA was extracted from 50 μL of whole blood collected in PAXgene tubes using custom Galenvs Total RNA kit (Galenvs) following manufacturer’s recommendations. Briefly, whole blood aliquot is mixed with 50 μL PBS and introduced onto the cartridge for automated protocol or processed manually for extraction in tubes. The mixture was first combined with 20 μL Proteinase K, and mixed. Lysis/binding buffer was then added to the solution and incubated at 55 °C for 10 min. For manual extraction in tubes, a DynaMag magnetic rack (#12321D Thermo Fisher Scientific) was used to capture magnetic nanoparticles. Following the capture of the RNA bound to the beads, two consecutive wash steps are performed. Elution was performed in 25 μL of nuclease-free water (Sigma-Aldrich). On-chip extraction of total RNA was performed using the automated protocol (Supplementary Fig. 5A) implemented on the centrifugal platform with the same reagents and volumes as for the manual extraction. For the on-chip capture of magnetic nanoparticles (MNPs), the external magnetic field was provided by a nickel-plated neodymium alloy disk magnet (D201, 1/8” in diameter, 1/32” in thickness; K&J Magnetics) which remained inserted in the designated area on the cartridge for the entire duration of the automated assay. The extracted RNA was assessed using NanoDrop One spectrophotometer and A260/A280 values were between 1.75 and 2.25, with typical yields in the range of 6–7 ng/µL for both manual and automated protocols. The extracted RNA was subsequently used in downstream RT-qPCR for assessment of RNA extraction efficiency as well as in on-chip ddPCR for determination of transcript copy number. Ct values for the 8-gene signature were similar for manual and automated extractions (ΔCt <1) with standard deviation slightly lower for the automated protocol.
qPCR
cDNA obtained from different patients were analyzed in a multiplexed qPCR using primer-probe sequences for genes of interest and housekeeping genes as internal controls for normalization. Each qPCR reaction consisted of 5 μL 10X PCR Buffer, 8 μL HotStar Taq Plus DNA Polymerase (Qiagen # 201205), 3 μL 25 mM MgCl2, 1 μL dNTPs, 5 μL 10X primer-probe mix (final concentration of 1 μM and 0.5 μM, respectively), 2 μL template, and 26 μL nuclease-free water (Sigma-Aldrich #W4502), for a total volume of 50 μL. Samples were tested in duplicate. A no-template control (NTC) reaction was included to assess for contamination. Thermal cycling was performed according to the manufacturer’s recommended protocol in a Bio-Rad CFX96 Touch Real-Time PCR Detection System (Bio-Rad). To quantify the copies of genes of interest each qPCR run included serial dilutions of cDNA (Takara) generating as such a standard curve. Cq values were plotted against the log concentration and linear regression was used to determine standard curves. The efficiency of each assay was 100 ± 10% and the R2 of each standard curve was >0.98.
RT-ddPCR
The ddPCR reaction master mix comprised of 5 μL 10X PCR Buffer (#201205, Qiagen), 8 μL HotStar Taq Plus DNA Polymerase, 8 μL 100X QuantiTect Virus RT Mix (#211015, Qiagen), 3 μL 25 mM MgCl2, 1 μL dNTPs, 5 μL 10X primer-probe mix (final concentration of 1 μM and 0.5 μM respectively), 2 μL template, and 18 μL nuclease-free water (Sigma-Aldrich), for a total volume of 50 μL. Template input was 2 µL of cDNA, RNA, or nuclease-free water for NTC samples. On-chip ddPCR assay was performed using an automated protocol implemented on the centrifugal platform (Supplementary Fig. 5B). Briefly, droplets containing template input in ddPCR reaction master mix were generated on-chip in fluorinated carrier oil (5% 00-8 FluoroSurfactant in HFE7500) (RAN Biotechnologies #008-FluoroSurfactant-5wtH-20G). The resultant emulsion was then transferred to the platform heater and cycled following manufacturer’s recommended protocol (20 min at 50 °C, followed by 5 min at 95 °C and 40 cycles of 15 s at 95 °C and 45 s at 60 °C, with ramp rate of 1 °C/s). Following thermal cycling, the emulsion was transferred to the chip for fluorescence imaging and data analysis. All experiments were performed in duplicate (no significant differences).
Microfluidic implementation of the SepsetER classifier detection process
The automated RNA extraction protocol (Supplementary Fig. 5A) starts with introduction of the sample in the RNA extraction chamber and installation of the cartridge on the platform. The software then executes a pre-programmed protocol sequence by initiating the platform to rotate. The first step of the automated workflow includes the transfer of a Proteinase K solution to the RNA extraction chamber, and bubble mixing. The lysis/binding buffer containing magnetic nanoparticles is subsequently transferred to the sample, mixed, and incubated for 10 min at 55 °C. The rotation speed is then increased to capture MNPs, and the lysate is transferred to the waste chamber. Two wash steps are then carried out sequentially by transferring the wash solutions from their respective chambers to the RNA extraction chamber. Finally, the purified RNA is eluted in the clean elution buffer.
To begin the cDNA synthesis and ddPCR protocol, a 2 µl aliquot of the eluted RNA is introduced on the ddPCR cartridge in the PCR mix chamber containing the RT-ddPCR master mix. Two cartridges, each having capacity to perform a single fourplex ddPCR reaction are operated in parallel to detect the 8-gene classifier. The automated sequence (illustrated in Supplementary Fig. 5B) commences by transferring the fluorinated oil into the droplet imaging chamber, followed by emulsification of RT-ddPCR master mix in the droplet generation chamber. The latter is performed by applying a positive pressure onto the ports of the master mix chamber to push the liquid through the resistive serpentine channel entering the array of nozzles connected to the shallow terrace merging into a deep reservoir of the droplet generation chamber. Upon completion of the droplet generation process, the rotation speed is reduced and positive pressure is applied in order to gently transfer the emulsion off-chip into the PCR tube located on the platform heater, using the world-to-chip interface. Following thermal cycling, the emulsion is transferred back on chip by applying positive pressure. The droplets sitting on the top of the oil in the neck of the imaging chamber are subsequently arranged in a monolayer suitable for imaging by applying a low negative pressure at ports of the oil reservoir. This step withdraws (back) the fluorinated oil from the imaging chamber into the oil reservoir and gently lowers the droplets into the shallower portion of the chamber. The pressure is slowly decreased to 0 psi until the monolayer formation is complete and the rotor is stopped, thus allowing the acquisition of fluorescence images for subsequent analysis.
Reporting summary
Further information on research design is available in the Nature Portfolio Reporting Summary linked to this article.
