Detection of Fusion Genes to Determine Minimal Residual Disease in Leukemia Using Next-Generation Sequencing

BACKGROUND
Measuring minimal residual disease (MRD), the persistence of leukemic cells after treatment, is important for monitoring leukemia recurrence. The current methods for monitoring MRD are flow cytometry, to assess aberrant immune phenotypes, and digital droplet PCR (ddPCR), to target genetic aberrations such as single-nucleotide variants and gene fusions. We present the performance of an RNA-based next-generation sequencing (NGS) method for MRD gene fusion detection compared with ddPCR. This method may have advantages, including the capacity to analyze different genetic aberrations and patients in 1 experiment. In particular, detection at the RNA level may be highly sensitive if the genetic aberration is highly expressed.


METHODS
We designed a probe-based NGS panel targeting the breakpoints of 11 fusion genes previously identified in clinical patients and 2 fusion genes present in cell lines. Blocking probes were added to prevent nonspecific enrichment. Each patient RNA sample was diluted in background RNA, depleted for rRNA and globin mRNA, converted to cDNA, and prepared for sequencing. Unique sequence reads, identified by unique molecular identifiers, were aligned directly to reference transcripts. The same patient and cell-line samples were also analyzed with ddPCR for direct comparison.


RESULTS
Our NGS method reached a maximum sensitivity of 1 aberrant cell in 10 000 cells and was mostly within a factor of 10 compared with ddPCR.


CONCLUSIONS
Our detection limit was below the threshold of 1:1000 recommended by European Leukemia Net. Further optimizations are easy to implement and are expected to boost the sensitivity of our method to diagnostically obtained ddPCR thresholds.


Detection of Fusion Genes to Determine Minimal Residual Disease in Leukemia Using Next-Generation Sequencing
Eddy N. de Boer,* Lennart F. Johansson, Kim de Lange, Anneke G. Bosga-Brouwer, Eva van den Berg, Birgit Sikkema-Raddatz, and Cleo C. van Diemen BACKGROUND: Measuring minimal residual disease (MRD), the persistence of leukemic cells after treatment, is important for monitoring leukemia recurrence. The current methods for monitoring MRD are flow cytometry, to assess aberrant immune phenotypes, and digital droplet PCR (ddPCR), to target genetic aberrations such as single-nucleotide variants and gene fusions. We present the performance of an RNA-based nextgeneration sequencing (NGS) method for MRD gene fusion detection compared with ddPCR. This method may have advantages, including the capacity to analyze different genetic aberrations and patients in 1 experiment. In particular, detection at the RNA level may be highly sensitive if the genetic aberration is highly expressed.

METHODS:
We designed a probe-based NGS panel targeting the breakpoints of 11 fusion genes previously identified in clinical patients and 2 fusion genes present in cell lines. Blocking probes were added to prevent nonspecific enrichment. Each patient RNA sample was diluted in background RNA, depleted for rRNA and globin mRNA, converted to cDNA, and prepared for sequencing. Unique sequence reads, identified by unique molecular identifiers, were aligned directly to reference transcripts. The same patient and cell-line samples were also analyzed with ddPCR for direct comparison.
RESULTS: Our NGS method reached a maximum sensitivity of 1 aberrant cell in 10 000 cells and was mostly within a factor of 10 compared with ddPCR.

Introduction
Leukemia, a cancer of the blood and bone marrow, is a highly heterogeneous disease caused by somatic variants. In addition to aneuploidies, single-nucleotide variants (SNVs), and indels in NPM1, CEBPA, FLT3, KIT, TET2, and DNMT3A genes, among others (1, 2, 3), >100 chromosomal rearrangements have been identified that result in a large number of unique fusion genes (1). In leukemia patients, the current workflow for detection of these types of aberrations is a combination of karyotyping, fluorescence in situ hybridization, and molecular techniques including PCR-based methods, Sanger sequencing, and, increasingly, next-generation sequencing (NGS) (4). Knowing which genetic aberrations a patient carries is important for diagnosis; it is one element of the WHO risk classification system (together with cell morphology and cell surface markers) (5), and this information is used to monitor the effectiveness of treatment (1,3,6,7).
Disease recurrence is monitored at intervals throughout treatment by measuring minimal residual disease (MRD) (8)(9)(10)(11)(12)(13), defined as the measurable persistence of leukemic cells after treatment (14,15). Although it is possible to identify leukemic cells in bone marrow at a threshold of 5% using only morphology, a far lower threshold is needed to determine risk classification for relapse or remission and intervention (16). Consequently, MRD is usually measured using flow cytometry, digital droplet PCR (ddPCR), or a combination of both (9,10,13,17,18). Flow cytometry distinguishes abnormal leukemic blast populations from normal progenitors using aberrant immunophenotypic features on the cell surfaces (19,20). In contrast, ddPCR allows absolute quantification of specific RNA or DNA sequences using fractioning, amplification, and probe-based labeling of RNA or DNA in thousands of droplets (21). For MRD detection at a genetic level, the patient-specific disease-causing variant is targeted. A more recently developed way of measuring MRD uses NGS to detect patient-specific aberrations at an RNA or DNA level (13,16,(22)(23)(24). This technique may have some advantages over ddPCR, such as the capacity to multiplex different aberrations and patients in 1 experiment (21,24). In addition, detection at the RNA level is preferable in the case of highly expressed fusion genes. This expression will result in more RNA copies of the fusion transcript compared with just 1 fusion DNA copy per cell; therefore, the sensitivity of RNA-based methods may be higher (1,23,24).
We developed an RNA-based NGS method for MRD detection using probes to enrich patient-specific fusion genes. We compared the sensitivity of our method with that of ddPCR, the current standard for MRD detection at a genetic level, applied to the same patient materials.

STUDY DESIGN AND LABORATORY WORKFLOW
We designed an enrichment panel based on gene fusions found in several patients and cell lines. To mimic the MRD situation, the panel was tested on dilutions of RNA samples of patients carrying the fusions in a background of control RNA. Our method for MRD detection is depicted in Fig. 1. In short, RNA is converted to cDNA, followed by end conversion, adaptor ligation, size selection, and amplification. For enrichment, fusion genes are hybridized with biotin-labeled probes and pulled out with streptavidin beads. Blocking oligonucleotides that target the wild-type transcripts and adapters are included during hybridization to prevent non-specific enrichment. The enriched fragments are then amplified and sequenced. Molecular identifiers are added to the adapters for identification of unique reads during data analysis.
Optimization of the procedure was performed in 2 steps to reach optimal sensitivity and specificity. To optimize enrichment, we tested different enrichment probe concentrations and different ratios of enrichment probe to blocking oligonucleotides. Size selection was then fine-tuned, and the number of PCR cycles was minimized. Finally, we compared the results of our optimized NGS procedure with ddPCR carried out on the same patient materials.

SAMPLE SELECTION, CELL STORAGE, AND RNA ISOLATION FROM STORED CELLS
We included 11 patients whose primary diagnosis was leukemia and who carried a specific fusion gene that was detected using the TruSight RNA Pan-Cancer panel (Illumina). We also included the cell lines ME-1 and Kasumi-1 (Leibniz Institute), which carry the RUNX1-RUNX1T1 and CBFB-MYH11 fusion genes, respectively (Supplemental Table 1). To obtain background control RNA, we mixed residual blood samples from patients with various indications, excluding leukemia.
The study protocol was approved by the ethics committee of the University Medical Centre Groningen (METC 2014.051, 10-2-2014), and all patients gave informed consent.

SAMPLE DILUTIONS
MRD samples from clinical patients contain a high percentage of wild-type transcripts mixed with a low percentage of fusion gene transcripts. The percentage of fusion gene transcript depends on the percentage of aberrant cells, the type of material, and the genes involved in the fusion. We simulated this situation by preparing mixes of RNA from 11 bone marrow samples from leukemia patients and 2 cell lines with background control RNA originating from human blood cells, while correcting for the percentage of aberrant cells in each patient sample. The percentage of aberrant cells was determined using fluorescence in situ hybridization or karyotyping (Supplemental Table 1). To lower processing costs, the number of samples was reduced by pooling 3 or 4 patient or cell-line dilutions that did not have overlapping fusion partners.

MRNA-TO-CDNA CONVERSION
To enrich for mRNA copies, 1000 ng of total RNA per diluted sample were depleted of globin RNA and rRNA using the deplete and fragment RNA part of the TruSeq stranded total RNA library prep globin protocol (Illumina). The mRNA was eluted without fragmentation and concentrated to a volume of 10 mL using Agencourt RNAClean XP beads (Beckman Coulter). The depleted and concentrated RNA solution was quantified using the Qubit RNA assay kit (ThermoFisher). Next, 70 ng of depleted RNA, originating from about 500 ng of input RNA, was converted to cDNA using the NEBNext Ultra II RNA library prep kit (New England Biolabs [NEB]). The fragmentation time was 15 minutes, and the extension time of the first-strand cDNA synthesis was increased to 50 minutes to ensure conversion of the longer RNA fragments. Purification after second-strand synthesis was performed using AMPure XP Beads (Beckman Coulter).

NGS PROCEDURE: PROBE DESIGN
We determined the highest expressed fusion gene transcript based on PanCancer expression data for all patients (Supplemental Table 2). The highest expressed fusion gene transcript in the cell lines had been determined in previous experiments (data not shown). We collected 120-bp sequences of the highest expressed fusion gene transcripts of the 13 translocations in a FASTA file (Supplemental Table 3). Care was taken to keep the exon-exon junction in the middle of the sequence. We selected housekeeping genes for quantification purposes and added the sequences of the most complete transcripts of AAAS, C10orf88, C12orf57, COX11, and DDX27 to the FASTA file using ENSEMBL GRCh37. Biotin-labeled probes were designed by Integrated DNA Technologies (IDT) to cover the sequences in the FASTA file using the xGEN LockDown probe protocol. The housekeeping genes were tiled with 120-bp biotin-labeled probes over the full transcripts.
To generate blocking-probes for wild-type transcripts, we collected 120-bp sequences covering all wildtype transcripts for each translocation partner in a FASTA file using ENSEMBL GRCh37 (Supplemental Table 3). Each sequence contained 100 bp of the exon present in the lowest expressed fusion gene transcript and 20 bp of the exon present in the highest expressed fusion gene transcript (Supplemental Table 2). We used only 20 bp of the highest expressed fusion gene sequence to prevent interruption of the enrichment procedure. IDT synthesized a 4-nmol/L solution of equimolarly mixed blocking probes for the regions in the FASTA file. Probe-quality criteria are described in Supplemental Table 3.

NGS PROCEDURE: LIBRARY PREPARATION
Library preparation was performed using the NEBNext Ultra II RNA Library Prep Kit (NEB). NEB adapters were substituted for unique molecular identifier (UMI) TruSeq dual-index duplex adapters (15 mmol/L; IDT) to remove duplicate reads and to reduce the error rate during data analysis. USER enzyme steps were skipped. After adding adapters, size selection for an average insert size of 400 bp was performed using AMPure XP Beads (Beckman Coulter). Library amplification and amplicon measurement were performed as described previously (25).

NGS PROCEDURE: ENRICHMENT PROCEDURE AND OPTIMIZATION
Targeted enrichment of the mixed sample dilutions was performed following the manufacturer's instructions (IDT; xGEN LockDown probes). We used 2 different capturing-probe concentrations: 0.0074 pmol/mL (manufacturer-delivered concentration) and a concentration of 0.00097 pmol/mL, which was proven to be successful in a previous study using a similar approach (25). Blocking probes were added during hybridization to prevent wild-type capturing, and several capturing-toblocking oligonucleotide ratios were tested to balance the final probe mix. Capturing-probe concentrations of 0.00097 pmol/mL were tested without blocking probes and in a dilution series with blocking-probe concentrations of 0.0024 pmol/mL (2.5:1 ratio), 0.012 pmol/mL (12.5:1 ratio), and 0.060 pmol/mL (62.5:1 ratio). Capturing-probe concentrations of 0.0074 pmol/mL were tested without blocking probes and in combination with blocking-probe concentrations of 0.019 pmol/mL (2.5:1 ratio) and 0.130 pmol/mL (17.5:1 ratio). We then used the best-performing composition for our final experiments. Eight separately indexed libraries were enriched in 1 mix.

NGS PROCEDURE: SEQUENCING AND DATA ANALYSIS
The sequence procedure was performed on a MiSeq (v2; 2 Â 150) or NextSeq (v2.5; 2 Â 150) sequencer (Illumina), following the manufacturer's instructions. To obtain a high number of unique reads per sample, we performed a NextSeq run. FASTQ files for index reads in MiSeq reporter were generated according to the manufacturer's instructions. FASTQ files for index reads using NextSeq were generated using bcl2fastq conversion software (Illumina). FASTQ files from separate NextSeq lanes were merged before data analysis.
Demultiplexing was performed automatically by MiSeq Reporter, and bcl2fastq conversion was performed using the unique sample bar codes. Data analysis was performed as described previously (25),omitting the optional prealignment step. In short, (1) FASTQ files were unzipped; (2) UMI sequences were extracted, and reads were put in an unmapped BAM file; (3) unmapped BAM files were converted to FASTQ files; (4) reads were aligned to reference FASTA files of the (fusion) genes of interest; (5) unmapped BAM files were sorted, and UMI information from the unmapped BAM files were connected to the mapped BAM files; (6) mapped reads were grouped by UMIs to create consensus reads; (7) BAM reads were converted to FASTQ for consensus reads, which were mapped to reference files of interest; and (8) unmapped consensus BAM files were sorted to merge UMI info into the mapped consensus BAM files.
The Integrative Genomics Viewer 2.3.1 (Broad Institute) (26,27) was used to view alignments using bam, bambai, fasta, and fasta.fai as input. The visibility of 1 aligned read with the translocation sequence is, in principle, enough for a positive detection result because the possibility that this happens by chance is negligible. Samtools view-c-F260 was used to estimate the number of reads aligned to the housekeeping genes. The number of unique indexes in an index file was counted using a custom script. Samtools view was used to count nonspecific capture of wild-type transcripts and genomic DNA.

DDPCR: PROBE AND PRIMER DESIGN
Probes and primers were designed using the manufacturer's recommendations (Bio-Rad). The Primer3 website (Whitehead Institute) was used for primer design, with care taken to keep the junction in the middle of the probe. Probes targeting the fusion were labeled with fluorescein amidite dye, and the control probe targeting the AAAS gene for quantification purposes was labeled with hexachloro-fluorescein dye. See Supplemental Table 4 for primer and probe characteristics.

DDPCR PROCEDURE
The iScript reverse transcription supermix for RT-qPCR (Bio-Rad) was used to convert 1-mg input RNA to cDNA. After optimizing the annealing temperature to 59 C, the ddPCR reaction was performed using the ddPCR supermix for probes (no dUTP; Bio-Rad) and an input of 150 ng RNA (separated into 3 wells). Quantasoft Pro software (Bio-Rad) was used for data analysis. The ddPCR was loaded conservatively to decrease the probability that there would be 2 fusion copies in 1 droplet and thus to increase the reliability of the extrapolation.

COMPARISON OF THE SENSITIVITY OF THE NGS METHOD WITH DDPCR
After determining the limit of detection of our NGS method, this dilution and a 10-times lower dilution were tested with ddPCR for each fusion gene. The number of fusion transcripts detected by ddPCR was divided by the number of fusion transcripts detected by NGS and converted to a factor of sensitivity by taking the dilution factors into account. In addition, we extrapolated the sensitivity to a ddPCR input of 500 ng for fair comparison to the NGS method.

MAXIMIZING THE PERFORMANCE OF THE NGS PROCEDURE BY BLOCKING WILD-TYPE CAPTURING
To maximize the specificity and sensitivity of fusion gene transcript enrichment in our procedure, we tested different capturing-probe concentrations (0.00097 and 0.0074 pmol/mL) and several blocking-to-capturing probe ratios. The optimal wild-type blocking result (78% wild-type expression blocked)-and, therefore, the highest specificity-was obtained with a 62.5:1 blocking-to-capturing probe ratio using the capturingprobe concentration of 0.00097 pmol/mL (Table 1 and  Supplemental Table 5). However, the sensitivity for fusion gene transcript capturing was higher with a 12.5:1 blocking-to-capturing probe ratio. Because maximal sensitivity was our priority, we used the 12.5:1 ratio for our final experiments. Increasing the capturing-probe concentration to 0.0074 pmol/mL had no effect on the sensitivity and a negative effect on blocking wild-type capturing ( Table 1 and Supplemental Table 5).

INCREASING SENSITIVITY BY DECREASING THE NUMBER OF PCR CYCLES
For fusion gene transcript detection, we reached the highest sensitivity by increasing the number of UMIs per sample, by lowering the number of pre-and postcapture PCR cycles, and by increasing the sequence output (Supplemental Tables 6 and 7). We assume that the transcription of the fusion gene has a limited influence on the total RNA concentration in patients; therefore, the RNA ratios of the patient RNA mixed with reference background RNA can be translated to cell ratios. As such, we reached a maximum sensitivity of 1:10 000 cells for NUP98-PSIP1 and NUP98-SET fusion gene transcripts ( Table 2).

COMPARISON OF SENSITIVITY OF THE NGS METHOD AND DDPCR
We reached a highest sensitivity of 1:100 000 cells using ddPCR. For an input of 150 ng RNA, the sensitivity for fusion gene transcript detection using our NGS method was mostly within a factor 10 of the ddPCR sensitivity (Table 2). However, the sensitivity for RUNX1-RUNX1T1 and PAX5-AUTS2 fusion gene transcript detection was, respectively, 255 and 120 times lower using NGS compared with ddPCR ( Table 2).

QUANTIFICATION OF TRANSCRIPTS USING HOUSEKEEPING GENES
As expected, the variance of the aligned reads of housekeeping genes within each sequence run decreased as the number of unique reads per sample increased. The normalized index-corrected reads that were aligned to the housekeeping genes are depicted in Supplemental Table  8. A higher number of unique reads per sample is caused by either decreasing the number of PCR cycles or increasing the sequence output per sample. In the most optimal situation, the coefficient of variation deviated between 0.05 and 0.12 (based on 8 samples per run; Supplemental Table 8).

Discussion
We developed an RNA-based NGS method for MRD fusion gene detection. To our knowledge, this study is the first time that a probe-based method has been tested for this purpose. We compared the results of our method and ddPCR for the same patients. We reached a maximum sensitivity of 1 aberrant cell in 10 000 cells for our method and 1:100 000 for ddPCR. Despite the lower sensitivity, our NGS method has advantages such as allowing analysis of multiple samples in a single experiment without the need to optimize each individual probe or to have prior knowledge of each individual target. We expect that expanding the panel to include other fusion genes and exonic regions will not influence the performance of the probes already included.
The required sensitivity of MRD is related to the type of treatment and the collection time (10). Although the literature is inconsistent, required sensitivity of 1:100 000 cells is regularly mentioned (10,17,28). We did not reach this level, but our detection limit was below the threshold of 1:1000 cells recommended by European Leukemia Net (16). In addition, the sensitivity of current reported methods is sometimes overestimated (10). Sensitivity at an RNA level is strongly influenced by the expression level of the fusion gene transcript in the cells of interest, and expression levels are strongly dependent on the type of fusion gene (13,24) and the type of material (i.e., bone marrow or blood) (1). Because of these different factors, if one wants to make fair comparison between methods, it is important to use the same samples for both methods. In general, The ddPCR method using 150 ng of input material showed 2-to 10-times higher sensitivity than our NGS method. However, the sensitivity of ddPCR was much higher for the RUNX1-RUNX1T1 and PAX5-AUTS2 fusion transcripts. This result might have been caused by the globin RNA and rRNA depletion, during which these low-abundance transcripts may have been depleted. Further experiments are needed to confirm this finding.
Although the required sensitivity of MRD detection is still the subject of debate, we tried to maximize it by blocking the capture of wild-type transcripts and increasing the number of unique sequence reads by decreasing the number of PCR cycles. However, further improvement is possible by blocking the capture of residual DNA fragments in the same way that we blocked the wild-type transcripts. In addition, further reduction in the number of PCR cycles will increase the number a NGS was carried out using an input of 500 ng of RNA; ddPCR was performed using 150 ng. After determination of the NGS limit of detection, this limit and a 10-times lower percentage were tested using ddPCR. The number of fusions was transformed to a factor of sensitivity. This factor was then extrapolated to a ddPCR input of 500 ng RNA. b Extrapolated from sample with a higher fusion percentage.
of unique reads in the sequence library (Supplemental Table 6 and 7). This can be achieved by increasing the amount of input RNA or expanding the captured region by including more housekeeping genes, additional fusion gene transcripts, or gene transcripts with SNVs. An additional benefit of adding gene transcripts with SNVs would be the ability to carry out MRD testing of different types of aberrations using 1 design. We expect that the increase in unique reads due to prevention of preferential amplification is higher than the additional reads needed to reach a given coverage because of the expansion of the design. MRD detection of aneuploidies would require another approach, unless a fusion or SNV were present in the same leukemic cell. The suggested modifications would probably lead to an NGS method that reaches a sensitivity similar to that of the current diagnostically used ddPCR method and with the aforementioned advantages of NGS.
In clinical MRD samples, the number of aberrant cells is unknown, in contrast to experimental setups. Therefore, in ddPCR methods, the fusion gene expression is related to the expression of a control gene for quantification purposes. In our NGS method, we showed that fusion gene expression can be compared with the constant expression of different housekeeping genes in a similar way. Given the variability in fusion gene expression between different samples, it is not possible to calculate the percentage of aberrant cells at an RNA level using either method.
A primer-based NGS method was published previously for MRD detection at an RNA level (24). This method reached a maximum sensitivity of 1 aberrant cell in 100 000 cells using cell dilutions, which was generally comparable to the sensitivity they reached with ddPCR and quantitative PCR. We note, however, that those authors used 2.5 times more input material for NGS, which may have affected the sensitivity. In its reported setup, their PCR-based method seems to be more sensitive than our probe-based method. However, PCR-and probe-based NGS methods each have unique pros and cons. PCR-based methods combine sample prep and enrichment in 1 assay, which makes them cheaper, less laborious, and faster than probe-based methods. In contrast, probe-based methods can handle larger amounts of input material and are more suitable for studying larger regions of interest. Both factors lead to a higher number of unique reads, which means that probe-based methods can reach higher sensitivity. In addition, larger enrichment probes are not susceptible to false-negative results caused by unexpected variants in the target sequence, and information about the exact sequence is less important. This is extremely relevant for the current diagnostic procedure, which usually still consists of methods like karyotyping and fluorescence in situ hybridization followed by Sanger sequencing to determine the exact breakpoints. With our probe-based method and the inclusion of all known fusion gene transcripts, the final Sanger sequencing step will not be necessary; however, we need to test this further.
In conclusion, our NGS method using RNA as input material allows analysis of multiple samples in a single experiment, without the need to optimize each individual probe or to know each individual target. Our method allows users to combine SNV targets and fusion detection in 1 design for MRD detection. Although our method requires technical improvements to reach higher sensitivity and to reduce sequencing costs, both goals are likely to be met in the near future. Furthermore, in contrast to current methods, our capture approach will make the use of a separate method to determine the exact sequence of a fusion target obsolete (21,23).

Supplemental Material
Supplemental material is available at Clinical Chemistry online.