Compositions and methods for identification of gene fusions

Abstract

The present invention relates to compositions and methods for identification of gene fusions. More specifically, the present invention relates to compositions and methods for ligating single stranded polynucleotides of unknown sequence. The resulting sequencing library, with inserts derived from the set of PCR amplified adapter-ligated DNA molecules are then subjected to computational analysis that includes filtering and error correction steps for the purposes of identifying among this set, any target gene with atypical exon configuration, or target gene which has been fused with any other DNA sequence such as a partner gene in a fusion configuration. Such genomic structural rearrangements or fusion events are clinically actionable with targeted inhibitors.

Description

COMPOSITIONS AND METHODS FOR IDENTIFICATION OF GENE FUSIONSFIELD OF INVENTION

[0001] The present invention relates to compositions and methods for identification of gene fusions. More specifically, the present invention relates to compositions and methods for ligating single stranded polynucleotides of unknown sequence and a method for identifying from a sequencing library of such polynucleotides any target gene with atypical exon configuration, or target gene which has been fused with any other DNA sequence such as a partner gene in a fusion configuration.BACKGROUND OF THE INVENTION

[0002] Malignancies are often driven by mutations, referred to as pathogenic variants, in the genomic DNA of cells. The most common types of pathogenic variants are single nucleotide variants within the coding regions of genes, however many other types of mutations can drive oncogenesis. Approximately 5% of solid tumors are driven by structural rearrangement of chromosomes which disrupt the linear order of the chromosome resulting in fusion of normally non-adjacent regions within or between chromosomes. This can result in gene fusions at the DNA level that are transcribed into RNA transcripts and translated into a chimeric protein with oncogenic properties.

[0003] Gene fusions, or translocations, is a class of genomic structural variants, and are estimated to contribute to disease in 16% of cancers and are the primary drivers of 1% of cancers1. Gene fusion prevalence, driver gene identities, and partner gene repertoires vary extensively between cancer types creating a need for partneragnostic pan-cancer gene fusion testing1-3. Gene fusion testing enables the use of targeted therapeutics available for treatment of kinase fusions4’ Multi-gene partneragnostic fusion gene detection tests that have high sensitivity and specificity at low cost are essential for driving uptake of gene fusion testing in countries with high cancer rates and low reimbursement5.

[0004] Such oncogenic fusions can involve fusion of a “driver gene” which provides oncogenic function with any of a large assortment of “partner genes”. A variable number of exons from the driver gene (often constituting a functional domain) is linked with a variable number of exons from the partner gene (zero exons are sometimes observed if only the promoter is required). The result of this fusion can bea chimeric protein with oncogenic function due to, for example, abnormal expression levels, abnormal cellular localization, or abnormal activity.

[0005] Gene fusions involving one or two genes commonly involve two genes with a truncated driver gene, often a kinase, that is fused with one of many partner gene promoters6. The orientation of gene fusions can vary with the driver gene at the 5’ or 3’ position, driving oncogenesis through abnormal expression levels, abnormal cellular localization, or abnormal activity, ultimately affecting pathways that affect tumor growth1 8. In the case of fusions involving a kinase, gene breakpoints can be located within intron regions that link a functional gene promoter with a driver kinase domain to produce a functional protein9. Other clinically significant structural changes include the loss of one or more exons within a single gene e.g. EGFR and MET exon 14 skipping events.

[0006] Molecular testing results guide the clinical management of cancer patients and the usage of targeted therapeutics available for treatment of gene fusion driven cancers. While gene fusions involve many functional classes of proteins, small molecule therapy development for oncogenic fusion genes has focused on kinase fusions Tyrosine kinase inhibitors (TKIs) have been developed over the last two decades as treatments for 5’ (FGFR and EGFR) and 3’ (ALK, CSF1, MET, NTRK, RET, ROST) receptor tyrosine kinase (RTK) gene fusions12 13. Gene fusions involving the ligand NRG1 also result in constitutive activity of ERBB family RTKs14. Similarly, MEK inhibitors show promise for BRAF and RAF1 gene fusions encoding serine / threonine receptor kinases15 16. The combination of fusion detection and the availability of tissue specific inhibitors has been transformative for patient outcomes with gene fusion driven cancers including prostate (ERG), bladder (FGFR3), lung (ALK, ROS1, NRG1), thyroid (RET), and bile duct (FGFR2)17. A new wave of tissue agnostic drugs such as those targeting NTRK gene fusions offers further hope for patients who previously required off-label therapeutic treatment13 17 18.

[0007] Current gene fusion detection methods include next generation sequencing (NGS) methods that can interrogate many targets for aberrations in combination with any partner gene exons, while whole transcriptome sequencing can interrogate the whole exome for gene fusions but can identify many chimeric reads or passenger gene fusions1 2. Therefore, classification methods can interrogate features that differentiate therapeutically relevant fusions from artifacts including driver gene identity and orientation, kinase domain retention, in-frame coding regions, and complete retention of fused exons30.

[0008] The detection of gene fusions is routinely performed on formalin fixed and paraffin embedded (FFPE) tissue using single gene histochemical methods including immunohistochemistry (IHC) and Fluorescence in-situ Hybridization (FISH)19. FFPE RNA may be fragmented which prevents ribo-depletion, or polyA tail based methods, and conversion into cDNA and amplification in general is inhibited by formalin- induced RNA modifications31~33. Additionally, the low fraction of degraded fusion transcripts in total RNA and the limited fraction of transcripts that are successfully extracted and converted into cDNA can be a challenge2434. A hybrid-capture approach for low quality FFPE RNA samples has been described but requires higher input template amounts35. FISH provides an indication of large-scale translocations that may result in a gene fusion but does not provide information about gene fusion identity or expression. IHC enables detection of aberrant expression of protein variants but can be challenging to interpret. Neither provides the base-level resolution produced by sequencing methods that enables pan-fusion detection and functional interpretation. Consequently, many clinical labs now utilize multi-gene panel amplicon and hybrid-capture methods with next generation sequencing (NGS) for improved sensitivity or lower costs and higher throughput 6,20,21 Both amplicon and hybrid capture assay technologies have unique challenges, however it is important for clinical assays cover the combinatorial diversity of partner and driver genes to guide clinical intervention3 1322. Therefore, assays that incorporate partner gene agnostic detection methods to detect any driver gene fusion are useful for clinical testing23~25, for example, for targets like ALK and NTRK that are promiscuous in terms of 5’ partners, have large introns, and are unpredictable in terms of recurrent breakpoints18. RNA-based panels enable partner-agnostic pan-fusion assessment due to the absence of introns which greatly reduces probe tiling, sequencing depth, or multiplexing requirements in addition to providing functional evidence of transcription23.

[0009] Polymerase chain reaction (PCR)-based methods have been used to exponentially amplify variants for identification or characterization, however such methods require some knowledge of the variant location and adjacent sequence to design flanking sequences that can be used as primers for targeted amplification. In the case of mutations involving structural rearrangements, the break site location is only loosely constrained in sequence space within both the partner and the driver genes resulting in a large potential sequence search space.

[0010] Partner-agnostic Reverse-Transcriptase (RT)-PCR-based methods provide an elegant strategy for fusion gene identification. First, cDNA (generated from total RNA) is used as a template in order to exclude large intronic regions thereby reducing sequence search space to only the shorter exons. Second, single primer extension from driver gene exons across exon junctions followed by amplification from a ligated universal adapter sequence removes the need for knowledge of partner gene sequence or identity. As a result, the sequence search space can be drastically reduced to only that of a select number of targeted driver gene exons. Subsequent PCR-amplification and sequencing of amplicons permits identification of atypical driver gene exon order or structural rearrangements involving fusion of driver gene exons with any other partner gene sequence.

[0011] Partner-agnostic PCR-based methods use a TA ligation method to add such an adapter to the ends of all DNA fragments from a patient sample to permit amplification from within the driver-gene region. Because TA ligation relies on a single nucleotide overhang this strategy has limited efficiency and conversion rates. As a result, sensitivity is decreased, and fusions that are poorly expressed, poorly represented, or have a limited number of amplifiable molecules can be difficult to identify. An alternative approach to improve ligation efficiency is to use longer overhangs that stabilize the interaction between ligation donor and acceptor molecules. This however requires a priori knowledge of the donor and acceptor sequences and limits ligation to between a single acceptor and donor molecule which makes this approach impractical for adapter ligation to a mixed population of fragments with unknown sequence as is the case in patient-derived cDNA fragments. It is also useful to analyse and process sequencing data from the PCR amplified adapter-ligated DNA molecules to identify specific genomic rearrangements or fusion events which may be clinically actionable. Therefore, a computer implemented method is useful to process the voluminous quantity of data so that results can be obtained in time to assist the patient.SUMMARY OF THE INVENTION

[0012] The present invention relates to compositions and methods for identification of gene fusions. More specifically, the present invention relates to compositions and methods for ligating single stranded polynucleotides of unknown sequence. The resulting sequencing library, with inserts derived from the set of PCR amplified adapter-ligated DNA molecules are then subjected to computational analysis.

[0013] In one aspect, the present invention provides a method for ligating single stranded DNA molecules by: i) providing a first single stranded DNA molecule comprising a known sequence of at least 10 nucleotides in length; ii) providing a second single stranded DNA molecule comprising an unknown sequence of least 10 nucleotides in length; iii) providing a single stranded RNA molecule comprising a 5’ end substantially complementary to the known sequence of the first single stranded DNA molecule and a 3’ end comprising a random sequence of least 10 ribonucleotides in length; iv) providing a DNA ligase; and v) combining the first single stranded DNA molecule, the second single stranded DNA molecule, the single stranded RNA molecule and the DNA ligase, where the first single stranded DNA molecule and the second single stranded DNA molecule are ligated to provide a ligated single stranded nucleic acid molecule when the unknown sequence of the first single stranded DNA molecule is substantially complementary to the random sequence of the single stranded RNA molecule.

[0014] In some embodiments, the first single stranded DNA molecule is blocked at the 3’ end.

[0015] In some embodiments, the first single stranded DNA molecule further includes a spacer sequence adjacent to the 3’ end of the known sequence.

[0016] In some embodiments, the known sequence of the first single stranded DNA molecule is from ten to about 100 nucleotides in length.

[0017] In some embodiments, the unknown sequence of the second single stranded DNA molecule is from ten to about 10,000 nucleotides in length.

[0018] In some embodiments, the single stranded RNA splint molecule is about 20 nucleotides in length and the 10 nucleotides at the 5’ end are substantially complementary to the known sequence of the first single stranded DNA molecule and the random sequence is substantially complementary to the unknown sequence of second single stranded DNA molecule.

[0019] In some embodiments, the first single stranded DNA molecule, the second single stranded DNA molecule, and the single stranded RNA molecule are combined under conditions suitable for annealing of the first single stranded DNA molecule and the second single stranded DNA molecule to the single stranded RNA molecule.

[0020] In some embodiments, the method further comprises use of an antiprimer molecule, where the first single stranded DNA molecule is combined with the single stranded RNA molecule and the second single stranded DNA molecule is combined with the antiprimer single stranded DNA molecule prior to combination of the first single stranded DNA molecule and the single stranded RNA molecule duplex and the second single stranded DNA molecule and the antiprimer single stranded DNA molecule duplex.

[0021] In some embodiments, the DNA ligase is a Chlorella virus DNA ligase.

[0022] In some embodiments, the second single stranded DNA molecule is in vitro ssDNA or single stranded cDNA.

[0023] In some embodiments, the second single stranded DNA molecule further includes a second known sequence 5’ to the unknown sequence.

[0024] In some embodiments, the second single stranded DNA molecule further includes a primer sequence at the 5’ end.

[0025] In some embodiments, the method further comprises amplifying the ligated single stranded nucleic acid molecule.

[0026] In some embodiments, the amplifying is performed by PCR.

[0027] In some embodiments, the PCR is linear PCR or nested PCR.

[0028] In an alternative aspect, the present invention provides a kit including: i) a first single stranded DNA molecule including a known sequence of at least 10 nucleotides in length; ii) optionally, a second single stranded DNA molecule including an unknown sequence of least 10 nucleotides in length; andiii) a single stranded RNA splint molecule including a 5’ end substantially complementary to the known sequence of the first single stranded DNA molecule and a 3’ end including a random sequence of least 10 ribonucleotides in length.

[0029] In some embodiments, the kit further includes a DNA ligase.

[0030] In some embodiments, the DNA ligase is a splint ligase.

[0031] In some embodiments, the kit further includes instructions for use.

[0032] In an alternative aspect, the present invention provides a computer implemented method for identifying from a plurality of read pairs R1 and R2 of a sequenced sample with inserts derived from a set of PCR amplified adapter ligated DNA molecules, any target gene fused with any other DNA sequence which are potentially clinically actionable, including: generating a k-mer index of predetermined length k from a reference genome or a transcriptome; mapping each sequencing read from the plurality of read pairs to exons of the reference genome with the k-mer index to identify corresponding exon regions of the sequencing read, and generating for each sequencing read a mapping vector including one or more blocks, where each block is either an exon block defined by a range of consecutive k-mer matches and their corresponding locations on an exon of the reference genome, or a mismatch block defined by a length of consecutive k-mer mis-matches; identifying exon blocks with overlapping consecutive k-mer matches in each of the mapping vectors, and merging the same exon blocks with adjustment of mapping coordinates or retaining the exon block having predetermined high priority criteria to resolve mapping ambiguities; combining both the mapping vectors corresponding to each read pair and sorting by genomic coordinates to provide intermediate mapping vectors, and merging any exon blocks with the same respective exons with adjustment of mapping coordinates to generate merged mapping vectors; identifying fusion junction candidates in the merged mapping vectors by searching for a pair of consecutive exon blocks from different genes, called fusionpartner genes, with a mismatch block in between, and recording the corresponding read pairs for each of the identified fusion junction candidates; clustering all the fusion junction candidates based on genomic coordinates of junction points of the pairs of consecutive exons from different genes where junctions closer than a predetermined number of base pairs are considered to be the same fusion junction candidate, and calculating a first predetermined set of metrics based on the read pairs including read depth; assembling a fusion transcript for each fusion junction candidate cluster using the corresponding recorded read pairs to generate assembled transcripts, and calculating a second predetermined set of metrics from the assembled transcripts; mapping the plurality of read pairs to the fusion transcripts, and calculating a third predetermined set of metrics based on the mapping output including a number of mapped reads; and categorizing any fusion transcript as being a high confidence reportable fusion when the orientation of both fusion partners is concordant, exactly one of the partners is a gene targeted by the assay, minimum thresholds for specific metrics of the first, second and third predetermined set of metrics are met, and the number of mapped reads and the read depth exceed predetermined thresholds.

[0033] In some embodiments, exon blocks having less than a minimum mapped bps are filtered out.

[0034] In some embodiments, the minimum mapped bps is 24 bps.

[0035] In some embodiments, each of the mapping vectors is represented by the expression (e_n, read_st— >read_end, map_st— >map_end) for each exon block and (miss, miss_st— >miss_end, len) for each mismatch block, where e_n is the exon unique identifier, read_st and miss_st are the starting base pair on the sequencing read, read_end and miss_end are the ending base pair on the sequencing read, map_st is the starting mapped region on the exon for read_st, map_end is the endingmapped region on the exon for read_end, miss is a status identifier for a mis-match, and len is the bps length of the mismatch block.

[0036] In some embodiments, the method further includes error correcting the mapping vectors to remove mismatch blocks interposed between two exon blocks with the same exon.

[0037] In some embodiments, error correcting includes removing the length of consecutive k-mer mis-matches of the mismatch block interposed between a first range of consecutive k-mer matches of a first exon block and a second range of consecutive k-mer matches of a second exon block of the same exon, when a continuation of the first range by the length of consecutive k-mer is sequentially followed by a beginning of the second range.

[0038] In some embodiments, error correcting to remove a mismatch block between two consecutive exon blocks with the same exon includes determining when read_end + len + 1 for the first exon block = read_st for the subsequent exon block to remove the mismatch block.

[0039] In some embodiments, error correcting to remove the length of consecutive k- mer mis-matches includes determining when map_end + len + 1 for a first exon = map_st for a subsequent exon to remove the length of consecutive k-mer mismatches.

[0040] In some embodiments, identifying fusion junction candidates includes executing an error correction applied on k-mer gaps g_1 , g_2 and g_k of the fusion junction candidates when g_k > 0, where g_1 is a distance between map_end of the first exon block to the corresponding exon extremity, g_2 is a distance between map_st of the second exon block to the corresponding exon extremity, g_k is the difference between the k-mer length k and the mismatch block len, and the k-mer gap correction steps include while g_k > 0 a. if g_1 > 0 theni. g_k = g_k - 1 ii. g_1 = g_1 - 1 b. else If g_2 > 0 then i. g_k = g_k -1 ii. g_2 = g_2 - 1

[0041] In some embodiments, the predetermined high priority criteria include an exon block with a larger mapping range than the other, an exon block where its exon comes from an assay targeted gene, an exon block matching a list of preferred candidates, or when both exons are from the same gene, selecting the exon with highest priority, calculated from sorting the corresponding reference genome transcripts for each exon in descending order of length.

[0042] In some embodiments, where if none of the predetermined high priority criteria is met, then remove the exon block having the largest exon gap, where an exon gap is a distance between map_end of a first exon and map_st of a subsequent exon to their respective exon extremities.

[0043] In some embodiments, the predetermined number of base pairs to cluster junction candidates is 50 or less.

[0044] In some embodiments, the method further includes determining a unique amplicon end by retrieving the last exon block mapping coordinate map_end and adding to it the length of the last mismatch block for each of the R2 reads.

[0045] In some embodiments, identifying fusion junction candidates further includes identifying normal junction candidates in the merged mapping vector by searching for pairs of consecutive exon blocks from a same gene with a mismatch block in between, and recording the corresponding read pairs for each of the identified normal junction candidates.

[0046] In some embodiments, the first predetermined set of metrics includes a number of unique amplicon ends, an estimated length of the fusion transcript, and a ratio between fusion junction candidates and normal junction candidates.

[0047] In some embodiments, the estimated length of the fusion transcript is determined by the average mapped length of all the read pairs in the corresponding junction candidate.

[0048] In some embodiments, assembling the fusion transcript includes the following steps for each fusion junction candidate cluster a. constructing a de-Bruijn graph with the k-mers of all cluster reads, b. executing bubble correction and tip removal algorithms, to correct sequencing errors, c. determining a number of contigs by traversing the de-Bruijn graph starting with a node that has the k-mer of the gene-specific primer that captured the fusion, and d. determining a longest contig as being a true assembled fusion transcript.

[0049] In some embodiments, the second predetermined set of metrics includes an assembled fusion transcript length.

[0050] In some embodiments, the third predetermined set of metrics includes the number of unique transcript ends determined from the mapping output.

[0051] In some embodiments, categorizing further includes categorizing any assembled fusion transcript as being an artifact if either the orientation of the fusion partners is not concordant or both partners are genes not targeted by the assay.

[0052] In some embodiments, categorizing further includes categorizing any assembled fusion transcript as being filtered if minimum thresholds for any one of the specific metrics is not met.

[0053] In some embodiments, categorizing further includes categorizing any assembled fusion transcript as being a supporting fusion if another assembled fusion with more reads is found on the same partner genes.

[0054] In some embodiments, categorizing further includes calculating a fusion score for each assembled fusion transcript based on any number of metrics of the first, second and third predetermined set of metrics, and the scoring function is calculated 'n'li — l-i . tor i G metricsusing the expression , where m_ij is the value of the metric / for the assembled fusion transcript j, wj is the weight associated with metric / , and p_\ is the expected mean of metric / , calculated from a control dataset with known fusions, comparing the fusion score to a first predetermined threshold and a second predetermined threshold, and categorizing the assembled fusion transcript as medium confidence if the fusion score is greater than the first predetermined threshold, categorizing the assembled fusion transcript as low confidence if the fusion score is less than the first predetermined threshold but greater than the second predetermined threshold, or categorizing the assembled fusion transcript as filtered if the fusion score is less than the second predetermined threshold, and outputting a result where low and medium confidence are subjected to orthogonal testing.

[0055] In some embodiments, the metrics used for the scoring function include the read depth, the number of unique reads, the assembled fusion transcript length, and a ratio of read depth between fusion junction candidates and normal junction candidates with the same gene as the fusion target gene.

[0056] This summary of the invention does not necessarily describe all features of the invention.BRIEF DESCRIPTION OF THE DRAWINGS

[0057] These and other features of the invention will become more apparent from the following description in which reference is made to the appended drawings wherein:

[0058] FIGURES 1A-D show a schematic representation of an DNA:RNA duplex formed by the the ssRNA splint, ssDNA phosphodonor (adapter) and ssDNA phosphoacceptor (amplicon) molecules. The following elements are represented: defined ribonucleotides (5’ ssRNA splint), random synthesized ribonucleotides (3’ ssRNA splint), defined deoxyribonucleotides (ssDNA phosphodonor and ssDNA phosphoacceptor), undefined (unknown) deoxyribonucleotides (ssDNA phosphoacceptor), C3 blocked deoxyribonucleotide (ssDNA phosphodonor), phosphorylated deoxyribonucleotide (ssDNA phosphodonor), defined deoxyribonucleotides used as primer sequence (ssDNA phosphoacceptor), defined deoxyribonucleotide (ssDNA antiprimer). A. The annealed phosphodonor hybrid duplex consisting of ssRNA splint and ssDNA adapter. B. The annealed phosphoacceptor duplex consisting of ssDNA amplicon and ssDNA antiprimer. C. The annealed ligation complex consisting of adapter and amplicon coordinated by RNA splint. D. The final ligation product after ligation by DNA ligase.

[0059] FIGURES 2A-C show oligonucleotide sequences used in adapter ligation. A. Nucleotide sequence of the ssDNA phosphodonor oligonucleotide. B. Ribonucleotide sequence of the ssRNA splint oligonucleotide. C. The ssDNA phosphodonor oligonucleotide and ssRNA splint oligonucleotide in annealed configuration showing base pairing. Randomly incorporated nucleotides are shown as italic A / . ssDNA phosphodonor modified bases (5’ phosphorylation, 3’ C3 spacer) are shown as underlined nucleotides A and C

[0060] FIGURE 3 shows a schematic representation of sequencing read orientation with respect to produced sequencing library molecules. P5 and P7 represent added Illumina sequencing platforms. Grey represents a library insert including ligated phosphodonor and phosphoacceptor molecules. The phosphoacceptor sequence range includes a central unknown sequence that can be determined by sequencing and a flanking known sequence used to orientate the sequencing data.

[0061] FIGURE 4 is a flowchart showing the computer implemented method for identifying target genes fused with any other DNA sequence which are potentially clinically actionable, according to a present embodiment

[0062] FIGURES 5A-B show examples of k-mer indexing and alignment. A. k-mer indexing. B. k-mer alignment.

[0063] FIGURE 6 shows an example of error-correction after a pseudo-alignment step according to a present embodiment;

[0064] FIGURE 7 shows an example of merged mapping of an R1 and R2 read, according to a present embodiment;

[0065] FIGURE 8 shows a junction signature example with k=5;

[0066] FIGURE 9 shows examples of alterations on the fusion gap g in the presence of sequencing errors / SNVs;

[0067] FIGURE 10 is a flowchart showing the fusion confidence generation method, according to a present embodiment; and.

[0068] FIGURE 11 shows example scores for fusion calls made on a control dataset with known fusions.

[0069] FIGURES 12A-G show ligation induced shift of dsDNA geneblock derived linear amplicon. (A) Anti-primers; (B-G) Visualization of 1 ng of components as follows: (B) Linear ssDNA template (C) With ligation and antiprimer (D) With ligation without antiprimer (E) Without ligation (no adaptersplint) (F) Without ligation (no ligase) (G) dsDNA template. The dsDNA ladder does not accurately size migration of non-dsDNA molecules.

[0070] FIGURES 13A-E show schematic representations of a diverse population of transcripts that can be converted into Illumina NGS libraries via ligation and amplification. (A) dsDNA template is amplified by PCR with a single gene specific primer (GSP). (B) The GSP is extended 5’-> 3’ by linear PCR to generate a ssDNA amplicon and the randomly located amplicon 3’ end is ligated to an adapter using the ligation method. (C) Successfully ligated amplicon is used as template for PCR amplification with adapter-specific primer and a second set of internal gene specific primers. (D) The internal GSP and adapter specific primers have Nextera tails that are used to add Illumina NGS library P5 and P7 sequences by PCR. (E) The final NGS library includes the entire population of initial amplicons with variable lengths and unknown 3’ regions.

[0071] FIGURES 14A-D show that geneblock ligation product derived NGS library is amplified from by PCR. (A) Library yields from template negative ligations. (B) Library yields from template positive ligations. (C-D) Visualization of 1 ng of components with Agilent Bioanalyzer High Sensitivity DNA Chip electrophoresis aselectropherogram and gel image composite (C) Specific product in ligation positive samples at 450bp. (D) Non-specific products present in ligation negative samples due to lack of template causing non-specific primer-primer interactions.

[0072] FIGURES 15A-D show that RNA transcript derived NGS library can be amplified from ligation product by PCR. (A) library yields from RNA transcript derived linear amplicon. (B-D) Visualization of 1 ng of components with Agilent Bioanalyzer High Sensitivity DNA Chip electrophoresis as electropherogram and gel image composite. (B) Non-specific products from template negative ligation. (C) Nonspecific products from ligase negative ligation. (D) Specific product in ligation positive sample.

[0073] FIGURE 16 shows the diversity of TUBB gene amplicon template molecule lengths as determined by paired end read sequencing of library fragments. The template molecule length is defined as the total length of the TUBB transcript between the start of sequencing read 1 (R1) and the end of sequencing read 2 (R2). For each unique sequence length, the count of unique sequences of that length is also shown.

[0074] FIGURE 17 is a flowchart showing an overview of the assay protocol. A. RNA is extracted and converted into cDNA. B. The cDNA is enriched for target gene exons by linear multiplex PCR, adapter ligated by random splint mediated ligation, and amplified for target gene exons by nested multiplex PCR. C. Libraries are generated by PCR and paired-end sequenced. D. Fusion calling pipeline flowchart. Reads from demultiplexed sample fastQ data files are filtered for quality, indexed using a k-mer sliding window and aligned against a reference exome reference to identify on-target reads. Candidate gene fusions are identified, evaluated, and assigned a call status. Starting at the top of the flowchart in part D of FIG. 17, the box of the first row represents an input file provided to the computing system. The boxes in the second row, and of the fourth to seventh rows represent processing steps executed by the computing system, and the boxes in the third and eighth rows represent output files provided by the computing system.

[0075] FIGURE 18 shows the TUBB gene exons 3 and 4. Enrichment and amplification primers are situated within exon 4, proximal to the exon junction, to enrich and amplify the adjacent exon, which in the normal TUBB transcript is exon 3.

[0076] FIGURE 19 shows a schematic representation of targeted molecule adapter ligation. The enrichment (white) and adapter (black) molecules are coordinated in a duplex by an RNA splint molecule (grey) with a random 10 nt 3’ extension. The enrichment molecule is duplexed at the 5’ end with a blocking oligo to prevent off- target ligation.

[0077] FIGURE 20 shows a schematic representation of the final composition of the target amplified and indexed library molecule. The target amplification molecule is generated by amplification with target amplification and adapter primers with Nextera (Illumina) 5’-tails which are converted into UDI indexed sequencing platforms using Nextera indexing primers (Illumina). The final library molecule contains the amplification primer site, exon junction site, and adjacent unknown sequence representing the expected adjacent exon in normal reads or an unexpected exon in gene fusions.

[0078] FIGURE 21 shows the read coverage of the TUBB gene transcript. Targetspecific enrichment was directed by an enrichment primer. After adapter-ligation, targeted amplification was performed using tailed amplification and adapter primers. Library conversion uses Nextera tails which permits paired-end sequencing into the insert from both directions. Unique reads are established based on their adapter-read start site positions relative to the aligned transcript.

[0079] FIGURE 22 shows on-target unique reads from three control samples in five runs: No template control (NTC), fusion-negative sample, and 18 fusion-positive sample.

[0080] FIGURE 23 shows counts of total unique fusion reads for each of the 18 fusion-positive sample fusions from ten observations for each scale of total template input (12.5ng, 25ng, 50ng, 100ng) plotted against manufacturer provided ddPCR expected copies for each fusion.

[0081] FIGURE 24 shows on-target unique fusion reads (median) for ten observations of the 18 fusion-positive sample with increasing input mass.

[0082] FIGURE 25 shows on-target unique fusion supporting read counts relative to increasing template input mass for 5 FFPE cell lines, and 15 FFPE archival clinical specimens (7: <2 years old, 8: 2-5 years old).

[0083] FIGURE 26 shows on-target unique fusion supporting read counts relative to qPCR CT value for template input mass groups of 20 FFPE samples. Samples with no detectable amplification were assigned a qPCR value of 40.DETAILED DESCRIPTION

[0084] The present disclosure provides, in part, methods for ligating single stranded polynucleotides that may be of unknown sequence, and a method for identifying from a sequencing library of such polynucleotides any target gene with atypical exon configuration, or target gene which has been fused with any other DNA sequence such as a partner gene in a fusion configuration.

[0085] The methods described herein permit sequence overhangs of at least 10 nucleotides Accordingly, the length of the sequence overhangs employed in the methods described herein provide greater stability to the interaction between ligation donor and acceptor molecules.

[0086] In some embodiments, the methods described herein do not require a priori knowledge of acceptor sequences. Accordingly, the methods described herein may be used in applications involving unknown sequences, for example, in adapter ligation to a mixed population of fragments with unknown sequence such as patient- derived cDNA fragments.

[0087] In some embodiments, the methods described herein are not limited to ligation to between a single acceptor and donor molecule. Accordingly, the methods described herein may be used in applications where multiple sequences can be ligated in a single assay, for example for ligation of a large population of unknown DNA molecules using a large population of random RNA splint molecules.

[0088] In one embodiment, there is provided a method for ligating single stranded deoxyribonucleic (DNA) molecules including: i) providing a first single stranded DNA molecule including a known sequence of at least 10 nucleotides in length; ii) providing a second single stranded DNA molecule including an unknown sequence of least 10 nucleotides in length;iii) providing a single stranded RNA molecule including a 5’ sequence substantially complementary to the known sequence of the first single stranded DNA molecule and a 3’ end including a random sequence of least 10 ribonucleotides in length; iv) providing a DNA ligase; and v) combining the first single stranded DNA molecule, the second single stranded DNA molecule, the single stranded RNA molecule and the DNA ligase, where the first single stranded DNA molecule and the second single stranded DNA molecule are ligated to provide a ligated single stranded nucleic acid molecule when the unknown sequence of the first single stranded DNA molecule is substantially complementary to the random sequence of the single stranded RNA molecule.

[0089] As used herein, the first single stranded DNA molecule may be referred to as the “adapter,” “ligation phosphodonor,” “ligation donor,” “phosphodonor” or “donor” molecule, and include a phosphate group at the 5’ end. The known sequence of the first single stranded DNA molecule may be a single defined sequence, for example, a sequence that can include one or both of the following characteristics: the two bases at the 5’ terminus may independently be adenine or thymine (A or T); the sequence may have a guanine or cytosine (G,C) ratio of approximately 50%. In some embodiments, the known sequence of the first single stranded DNA molecule may be 5’-ATCTGTCTCTTATACACATCTCCGAGCCCACGAGAC-3’ (SEQ ID NO: 1). It is to be understood that multiple adapter sequences can be used in a single reaction, provided that the reaction includes single stranded RNA molecules including a sequence complementary to the sequence of each adapter.

[0090] In general, the known sequence of the first single stranded DNA molecule is of a length sufficient to allow annealing to the substantially complementary 5’ sequence of the single stranded RNA molecule to, for example, form a stable duplex. It is to be understood that amount of overlap between the molecules determines the stability and specificity of the ligation duplex, as known in the art. For example, shorter or longer overlaps affect the temperature at which the ligation needs to be performed as shorter overlaps require lower temperatures for stabilization of short duplexes while longer overlaps require higher temperatures to achieve acceptable specificity. It is also to be understood that melting temperatures for DNA:RNA hybrid duplexes are more stable for a given length of complementarity than DNA:DNA. By“substantially complementary” is meant a sequence that is capable of annealing / hybridization under conditions, for example, conditions as described herein or known in the art to form a stable duplex. A “stable duplex” is a molecule that has at least a 4 base pair overlap between the single stranded RNA molecule and each of the first and second single stranded DNA molecules. In some embodiments, a “stable duplex” has at least a 10 base pair overlap between the single stranded RNA molecule and each of the first and second single stranded DNA molecules i.e., for a total of 20 base pair overlap.

[0091] In some embodiments, the known sequence of the first single stranded DNA molecule is of a length suitable for use as a primer.

[0092] The known sequence of the first single stranded DNA molecule may be at least 10 nucleotides in length. In some embodiments, the known sequence of the first single stranded DNA molecule may be 10 nucleotides to about 100 nucleotides in length, or any value in between and inclusive of this range, for example, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70, 75. 80, 85, 90, 95 or 100 nucleotides in length.

[0093] In some embodiments, the known sequence of the first single stranded DNA molecule may be blocked at the 3’ end. A DNA sequence may be “blocked at the 3’ end” if it includes a modification that prevents 3’ extension, for example, during PCR (polymerase chain reaction). Examples of such 3’ modifications include, without limitation, phosphorylation, spacers, modified bases, phosphorothioate bonds, fluorophores, such as C3 Spacer phosphoramidite, Dideoxycytidine (ddC), or 3' Inverted dT.

[0094] In some embodiments, the known sequence of the first single stranded DNA molecule may include a spacer sequence adjacent to the 3’ end of the known sequence i.e., between the known sequence and the 3’ block. The spacer sequence may be a low GC sequence. By “low GC” sequence is meant a sequence with fewer than 50% of the constituent nucleobases being guanine or cytosine (G,C). In some embodiments, the low GC sequence may be two nucleobases, of which one or both may be adenine or thymine (A, T).

[0095] As used herein, the second single stranded DNA molecule may be referred to as the “template,” “ligation phosphoacceptor,” “ligation acceptor,” “phosphoacceptor” or “acceptor” molecule, and include a hydroxyl group at the 3’ end. The unknown sequence of the second single stranded DNA molecule may be obtained from anysuitable source or sample. For example, the unknown sequence of the second single stranded DNA molecule may be any in vitro ssDNA or single stranded cDNA. In some embodiments, the unknown sequence of the second single stranded DNA molecule may be obtained from a patient known to have, or suspected of having a cancer, for example, of blood, lung, breast, or colorectal cancer. In some embodiments, the unknown sequence of the second single stranded DNA molecule may be obtained from a patient known to have, or suspected of having a cancer with an oncogenic structural rearrangement in a target gene, such as a kinase gene fusion. In some embodiments, the gene fusion may be, without limitation, a 5’ (for example, FGFR or EGFR) or a 3’ (for example, ALK, CSF1 , MET, NTRK, RET, or ROS1) receptor tyrosine kinase (RTK) gene fusion; a BRAF or RAF1 gene fusion encoding a serine / threonine receptor kinase; a gene fusion involving the ligand NRG1 ; a gene fusion driven cancer including prostate (ERG), bladder (FGFR3), lung (ALK, ROS1 , NRG1), thyroid (RET), or bile duct (FGFR2); or a NTRK gene fusion. In some embodiments, the unknown sequence of the second single stranded DNA molecule may be obtained from a patient known to have, or suspected of having, a cancer that is a solid tumour. In some embodiments, the cancer may be one where there are therapeutic drugs available or in development and the methods of described herein may be used to determine appropriate therapies for the patient. In some embodiments, the therapeutic drug may be, without limitation, a tyrosine kinase inhibitor (TKI) for the treatment of 5’ (FGFR or EGFR) and 3’ (ALK, CSF1 , MET, NTRK, RET, ROS1) receptor tyrosine kinase (RTK) gene fusions, a MEK inhibitor for the treatment of BRAF or RAF1 gene fusions encoding serine / threonine receptor kinases, etc. In some embodiments, the patient may be a human.

[0096] As used herein, a “target gene” is a gene known to have, or suspected of having, a structural rearrangement of chromosomes which disrupt the linear order of the chromosome resulting in fusion of normally non-adjacent regions within or between chromosomes that can result in gene fusions at the DNA level that are transcribed into RNA transcripts and translated into a chimeric protein with oncogenic properties. Accordingly, a target gene may include an “oncogenic fusion” which involves fusion of a “driver gene” which provides oncogenic function with any of a large assortment of “partner genes”. In an oncogenic fusion, a variable number of exons from the driver gene (which may include a functional domain) can be linked with a variable number of exons from the partner gene or, alternatively, with no exons from the partner gene if, for example, only the promoter is required, resulting in a chimeric protein.

[0097] In some embodiments, the target gene may be, without limitation, kinase encoding genes involved in solid tumour oncogenesis or genes involved in fusions in blood cancers. In some embodiments, the target gene may be, without limitation, selected from the genes listed in Table 1.

[0098] In general, the unknown sequence of the second single stranded DNA molecule is of a length sufficient to allow annealing to the random 3’ sequence of the single stranded RNA molecule to, for example, form a stable duplex, as described herein or known in the art.

[0099] In some embodiments, the unknown sequence of the second single stranded DNA molecule may be at least 10 nucleotides in length. In some embodiments, the unknown sequence of the second single stranded DNA molecule may be 10 nucleotides to about 10,000 nucleotides in length, or any value in between and inclusive of this range, for example, 10, 50, 100, 150, 200, 250, 300, 350, 400, 450, 500, 600, 700, 800, 900, 1000, 1500, 2000, 2500, 3000, 3500, 4000, 4500, 5000, 5500, 6000, 6500, 7000, 7500, 8000, 8500, 9000, 9500, or 10000 nucleotides in length.

[0100] In some embodiments, the second single stranded DNA molecule may further include a second known sequence 5’ to the unknown sequence. The second known sequence of the second single stranded DNA molecule may be obtained from a target gene, as described herein.

[0101] In some embodiments, the second single stranded DNA molecule may include a primer sequence 5’ to the second known sequence. In some embodiments, the primer sequence may be provided at the 5’ end of the second known sequence. In some embodiments, the second known sequence may be the primer sequence. In some embodiments, the primer sequence may be obtained from a target gene, as described herein, for example, in Table 3 (Outer Primers).

[0102] In some embodiments, the second single stranded DNA molecule including the unknown sequence and the second known sequence, and the primer sequence may be referred to as an “amplicon.”

[0103] Single stranded DNA molecules may be synthesized in vitro ssDNA oligonucleotides, cDNA molecules generated by reverse transcription of RNA, or generated by linear amplification from double stranded DNA, such as cDNA. It is tobe understood that any suitable method for generating single stranded DNA may be used.

[0104] As used herein, the single stranded RNA molecule may be referred to as the “splint” molecule. In general, the 5’ sequence of the single stranded RNA molecule may be substantially complementary to the known sequence of the first single stranded DNA molecule and of a length sufficient to allow annealing to the known sequence of the first single stranded DNA molecule to, for example, form a stable duplex as described herein or known in the art.

[0105] In some embodiments, the 5’ sequence of the single stranded RNA molecule may be at least 10 nucleotides in length. In some embodiments, the 5’ sequence of the single stranded RNA molecule may be 10 nucleotides to about 100 nucleotides in length, or any value in between and inclusive of this range, for example, 10, 15, 20, 25, 30, 35, 40, 45, 50, 55, 60, 65, 70, 75. 80, 85, 90, 95 or 100 nucleotides in length.

[0106] By “random sequence” is meant a DNA sequence that contains the DNA nucleotides in no particular order. The random sequence of the single stranded RNA molecule may be generated by any suitable means such that all possible permutations of the four nucleotides (A, T, G, C) are enabled. In some embodiments, 4A10 different permutations of the 3’ random sequence region are contemplated.

[0107] In general, the random sequence of the single stranded RNA molecule may be of a length sufficient to allow annealing to the unknown sequence of the second single stranded DNA molecule to, for example, form a stable duplex as described herein or known in the art. In some embodiments, the random sequence of the single stranded RNA molecule may be at least 10 nucleotides in length. In some embodiments, the random sequence of the single stranded RNA molecule may be 10 nucleotides to about 20 nucleotides in length or any value in between and inclusive of this range, for example, 10, 11 , 12, 13, 14, 15, 16, 17, 18, 19, or 20 nucleotides in length.

[0108] In some embodiments, the single stranded RNA molecule may be about 20 nucleotides in length where 10 nucleotides at the 5’ end are substantially complementary to the known sequence of the first single stranded DNA molecule and 10 random nucleotides at the 3’ end are substantially complementary to the unknown sequence of the second single stranded DNA molecule.

[0109] In some embodiments, the sequence of the single stranded RNA molecule may be 5’- rGrUrCrUrCrGrUrGrGrGrCrUrCrGrGrArGrArUrGrUrGrUrArUrArArGrArGrArCrArGrAr U rNrNrNrNrNrNrNrNrNrN-3’ (SEQ ID NO: 2). It is to be understood that the sequence of the single stranded RNA molecule is unrestricted, as long as it has complementarity with the phosphodonor adapter sequence, for example, of at least ten nucleotides. The randomly incorporated ribonucleotides (‘rN’) are at the 3’ end and include at least 4 nucleotides.

[0110] The first and second single stranded DNA molecules may be mixed with the single stranded RNA molecule under conditions suitable for annealing of the first and second single stranded DNA molecules to the single stranded RNA molecule, as described herein or known in the art.

[0111] In some embodiments, the first single stranded DNA molecule may be mixed with the single stranded RNA molecule under conditions suitable for annealing of the first single stranded DNA molecule to the single stranded RNA molecule, and second single stranded DNA molecule may be mixed with an antiprimer molecule under conditions suitable for annealing of the second single stranded DNA molecule to the antiprimer molecule. The first single stranded DNA molecule / single stranded RNA molecule duplex may then be mixed with the second single stranded DNA molecule / antiprimer molecule under conditions suitable for annealing of the second single stranded DNA molecules to the single stranded RNA molecule. As used herein, “antiprimer” refers to a sequence that is the reverse complement of at least a portion of the second known sequence of the second single stranded DNA molecule. In general, the antiprimer molecule may be of a length sufficient to allow annealing to the second known sequence of the second single stranded DNA molecule to, for example, form a stable duplex as described herein or known in the art.

[0112] The first and second single stranded DNA molecules, when annealed to the single stranded RNA molecule, are positioned such that they are compatible for end-to-end ligation by a “DNA ligase,” which is defined herein as a ligase capable of such end-to-end ligation of two single stranded DNA molecules coordinated in positioned by a RNA splint molecule by way of hybrid duplex base-pairing, such as a PBCV-1 DNA Ligase (Chlorella virus DNA Ligase or SplintR® ligase). Ligation conditions may be determined based on, for example, the length and composition of the DNA and RNA molecules (e.g., the first and / or second single stranded DNAmolecules, the single stranded RNA molecule) that would provide sufficient annealing. Exemplary ligation conditions are described herein.

[0113] The ligated molecule can be amplified using a primer complementary to the adapter sequence, such as the Nextera i7-R (5’-GTCTCGTGGGCTCGGAG-3’; SEQ ID NO: 3) and a primer internal to the primer used in the amplicon sequence. The internal primer sequence may be obtained from a target gene, as described herein. The internal primer sequence may further include a common 5’ sequence, for example, the Nextera XT universal amplification sequence (5 - TCGTCGGCAGCGTCAGATGTGTATAAGAGACAG-3’ (SEQ ID NO: 4)). It is to be understood that any suitable platform can be used, for example, the Fluidigm CS1 / CS2 platform and associated sequences. In some embodiments, the internal primer sequences may be those set forth in Table 6.

[0114] Accordingly, in some embodiments, the first set of target specific primers used in linear amplification to generate the phosphoacceptor molecule is termed the “outer” primer set. Once the phosphoacceptor molecules are ligated to the phosphodonor adapter, then in combination with adapter-specific primer, an “inner” primer set is used to amplify target regions internal to the initial target-specific primer sites. For the initial linear amplicon, the first primer set used is on the outside, while the second primer set used is on the inside with respect to the molecule. The outer set is used for linear amplification, the inner set is used for ligation product PCR amplification. For example:Linear amplicon: 5’-‘OUTER GSP’-‘INNER GSP’-unknown-3’Ligation product: 5’-‘OUTER GSP’-‘INNER GSP’-unknown-adapter-3’PCR amplicon: 5’-‘INNER GSP’-unknown-adapter-3’

[0115] A ligation complex can be generated by annealing the products of two annealing processes, whereA is the adapter phosphodonor moleculeB is the RNA splint moleculeC is the amplicon phosphoaccepter moleculeD is the antiprimer moleculeAnnealing 1 : A + B - > ABAnnealing 2: C + D -> CDAnnealing 3 and ligation: AB + CD + ligase -> ABCD, where A and C are covalently linked by the ligase through the condensation reaction of a phosphodiester bond. B and D interact with the newly created AC molecule through non-covalent hydrogen bonding. Molecules B and D are degraded, lost, or removed, in subsequent manipulation by denaturation or purification and molecule AC is used for subsequent manipulation such as PCR.

[0116] Accordingly, the final ligation involves mixing an annealed hybrid duplex with an annealed DNA duplex, each of which have single stranded extensions (the RNA splint random sequence and the linear amplicon unknown sequence). The extensions of these two molecules can base pair to coordinate ligation by the DNA ligase, e.g., splint ligase. The splint molecule 3' ssRNA extension base pairs with 3’ ssDNA extension of amplicon to form hybrid duplex positioning 3’ terminal hydroxyl residue of amplicon DNA adjacent to 5’ terminal phosphate residue of adapter DNA. Before phosphodiester bond is formed these two ends are approximately centrally coordinated by splinting by the RNA molecule.

[0117] Antiprimer molecules may be used to block excess outer primer from the linear amplification step from participating in the ligation and to reduce the amount of non-specific artifact generated as a result of undesired ligation of excess primer. Accordingly, the addition of a limiting amount of antiprimer can result in base pairing and formation of a duplex between antiprimer and outer primer to form dsDNA duplexes which are not ligation templates. The linear amplicon is the result of linear extension of the outer primers. For these molecules, the 5’ sequence corresponding to the outer primer also participates in base pairing with the antiprimers, however the 3’ linear extension remains in ssDNA configuration that is available for ligation.

[0118] In some embodiments, amplification of the ligated molecule can be used to generate sequencing libraries, as described herein.

[0119] In some embodiments, the present disclosure provides a composition comprising:i) a first single stranded DNA molecule including a known sequence of at least 10 nucleotides in length; ii) a second single stranded DNA molecule including an unknown sequence of least 10 nucleotides in length; iii) a single stranded RNA molecule including a 5’ sequence substantially complementary to the known sequence of the first single stranded DNA molecule and a 3’ end including a random sequence of least 10 ribonucleotides in length.

[0120] In some embodiments, the present disclosure provides a composition comprising two or more of the sequences set out in each of Tables 3, 5 and / or 6. In some embodiments, the composition comprises ten or more of the sequences set out in each of Tables 3, 5 and / or 6. In some embodiments, the composition comprises twenty or more of the sequences set out in each of Tables 3, 5 and / or 6. In some embodiments, the composition comprises thirty or more of the sequences set out in each of Tables 3, 5 and / or 6. In some embodiments, the composition comprises forty or more of the sequences set out in each of Tables 3, 5 or 6. In some embodiments, the composition comprises all of the sequences set out in each of Tables 3, 5 and / or 6.

[0121] In some embodiments, the present disclosure provides a kit comprising the composition, as described herein, and a DNA ligase, such as a splint ligase, optionally with instructions for use.

[0122] In some embodiments, the present disclosure provides a kit comprising two or more of the primers set out in each of Tables set out in 3, 5 and / or 6, as described herein, optionally with instructions for use. In some embodiments, the kit comprises ten or more of the primers set out in each of Tables 3, 5 or 6. In some embodiments, the kit comprises twenty or more of the primers set out in each of Tables 3, 5 and / or 6. In some embodiments, the kit comprises thirty or more of the primers set out in each of Tables 3, 5 and / or 6. In some embodiments, the kit comprises forty or more of the primers set out in each of Tables 3, 5 and / or 6. In some embodiments, the kit comprises all of the primers set out in each of Tables 3, 5 and / or 6.

[0123] The ligated single stranded polynucleotide (including a primer sequence) can be used in amplification methods such as PCR, for example linear PCR and / or nested PCR.

[0124] The amplified products may be analysed using computational methods, as generally described in Figure 4, which sets out steps in the fusions calling pipeline.

[0125] In some embodiments, the methods described herein may be useful in analysing a sample containing poor quality RNA for example, a Formalin-Fixed Paraffin-Embeeded (FFPE) tissue specimen obtained from, for example, a clinical patient. In some embodiments, the methods described herein may be useful in analysing a sample with low cellularity, for example, a small biopsy or fine needle aspirate, that yield low amounts of RNA. In some embodiments, the methods described herein may be useful in analysing a sample with low tumor purity, for example, a sample with very low fusion to normal read ratios, as described herein. In some embodiments, in low tumor purity cases, macrodissection of the FFPE section after evaluation by a pathologist or deeper sequencing may be performed. In some embodiments, in addition to qPCR CT values, the methods described herein provide multiple metrics that enable contextualization of low confidence or negative fusion calls including, for example, on-target reads, unique transcript reads, total and unique fusion reads, and fusion-to-normal read ratios.

[0126] The present invention will be further illustrated in the following examples.

[0127] EXAMPLES

[0128] Example i

[0129] ssDNA ligation template phosphoacceptor oligonucleotides were generated by linear amplification from double stranded cDNA.

[0130] This mixed population of ligation acceptor molecules with unknown 3’ sequence is then mixed with an adapter molecule with a single defined sequence (ligation donor molecule). The donor and acceptor molecules are then coordinated end-to-end in a ligation compatible orientation with a splinting RNA oligo through formation of a DNA:RNA duplex depending on sequence-based complementarity (Figure 1).

[0131] A synthesized mixture of RNA oligos with sequence matching the reverse complement of the adaptor molecule followed by a 10nt undefined stretch of nucleotides with 4A10 sequence permutations is used. The RNA oligo therefore formsan RNA:DNA duplex with the adaptor leaving a 10nt ssRNA extension for interaction with a sequence compatible ligation acceptor molecule. This design permits any ssDNA molecule with any 10nt 3’ sequence to find a RNA oligo with compatible sequence allowing for formation of a DNA:RNA duplex and ultimately ligation of the phosphoacceptor ssDNA and phosphodonor ssDNA adapter molecule ends which are coordinated in to position by the ssRNA splint for formation of the phosphodiester bond by the DNA ligase.

[0132] The adapter ligated ssDNA molecules can then be amplified by traditional PCR, converted into a library for sequencing, sequenced, and interrogated for sequence and identification of any fusion genes (Figure 2).

[0133] Methods:

[0134] Gene and Exon Selection: Genes were selected to be included in targeted content if they contribute to established oncogenic structural rearrangements in lung, breast, or colorectal solid tumours and if there are available therapeutic drugs available or in development. Exons for these genes were selected to be included in targeted content if they had more than one reported incidence in the literature or in large public databases (for example, The Cancer Genome Atlas Program (TCGA), and the Catalogue of Somatic Mutations in Cancer (COSMIC), etc.).

[0135] Primer and Oligo Design Three sets of multiplexed oligos were used in this assay design, in addition to a random RNA splint pool, a universal adapter DNA oligo, and a DNA oligo for universal adapter amplification. The universal adapter and amplification oligo sequences employ Nextera XT amplification sequences for compatibility with existing library amplification methods. The RNA splint was designed to be the reverse complement of the adapter oligo with the addition of a low GC spacer sequence and a random 10mer sequence at the 3’ end. For the multiplexed primer pools, in general, the inner primers were designed to be as close as possible to 20bp downstream of the exon junction, and the paired outer primer designed to be as close as possible to the inner primer without overlapping. The anti-primer set represented the reverse complement sequences of all the outer primers. After design of the inner primers a Nextera XT universal amplification sequence was appended to the 5’ end of the oligo to provide a universal sequence platform for library amplification steps.

[0136] More specifically, targeted exons were amplified using two sets of primers - an “outer” primer distal to the expected fusion junction and an “inner” primer proximal to the expected fusion junction. All the primers for each set were combined to form an “outer” multiplexed pool and an “inner” multiplexed pool.

[0137] The “outer” multiplexed pool was used for targeted enrichment of cDNA template by primer-directed linear amplification (there is no primer partner primer for the alternate strand and therefore amplification is linear rather than exponential). After ligation of an adapter to the target-enriched cDNA, the “inner” multiplexed pool was used with an adapter-sequence specific general amplification primer for nested amplification of adapter-ligated target-enriched cDNA. Each ligation product involving a target transcript was amplified from one direction by the target primer and from the other direction by a common primer with sequence complementarity to the adapter added to target amplicon by ligation.

[0138] While the “outer” primer included only the gene-specific sequence, the “inner” primers included a 3’ gene-specific sequence in addition to a common 5’ tail which adds the Nextera XT universal amplification sequence for universal amplification. The target amplicon ligation product therefore possesses at each flanking end, universal amplification sequences (one added by the adapter and one added by the “inner” primer) which were used for amplicon indexing and addition of Illumina sequencing p5 / p7 platforms by PCR using the IDT for Illumina Nextera Unique Dual indexing primers.

[0139] To design the multiplexed oligos a custom software was developed which performs the following steps: 1 . Gene, transcript, and exon annotations are imported from Ensembl database via Biomart. 2. Single Nucleotide Polymorphisms (SNPs) are imported from dbSNP. 3. An existing primer panel listing exon targets and any primers, if already designed, is imported. 4. Primers are designed for all desired exon targets using primer3 with parameters including exon identity, panel context file, proximity to exon junction, GC content, oligo length, SNP filter, and undesired characteristic filter. Parameters are progressively relaxed to achieve a minimum number of candidate primers. 5. From all the candidate primers, inner and outer pools are assembled with a semi-supervised method that minimizes undesired primer-primer interactions by way of calculating all possible Gibb’s free energy values for possible pools and selecting the pool that minimizes the total sum of these values.

[0140] RNA Extraction: RNA was extracted from FFPE samples using a Maxwell RSC instrument (Promega) and RNA FFPE kit (Promega) according to manufacturer instructions with elution in 40uL EB.

[0141] RNA Quality Assessment: RNA was assessed using the Bioanalyzer RNA Pico assay (Agilent) and a DV200 ratio representing the relative fraction of molecules > 200nt was calculated.

[0142] cDNA Synthesis: First strand cDNA was generated from total extracted RNA with random hexamers using the Superscript IV system (Thermofisher) according to manufacturer instructions. Second strand cDNA synthesis was generated using the NEB second strand synthesis module according to manufacturer instructions. cDNA synthesis product was purified with a 1 ,8X Ampure XP SPRI purification and elution in 25uL EB.

[0143] cDNA Quality Assessment: High quality RNA was converted with 100% efficiency to DNA resulting in a 1 :1 conversion efficiency which is monitored. As a result of modification of RNA molecules and the presence of inhibitors, RNA derived from FFPE material frequently has a conversion below 10%. Converted cDNA was assessed for amplifiability with a real-time PCR based assay according to manufacturer instructions with in-house validated threshold values (0.1) and cycle threshold (CT) values.

[0144] Target Enrichment: The assay employed a panel consisting of 19 target genes and 2 control genes (Table 1) covering 150 exons (Table 2).Table 1 : Target Genes

[0145] Table 2: Target Gene Exons

[0146] The 19 target genes listed in Table 1 represent oncogenic driver genes frequently observed in solid tumour structural rearrangements that have FDA approved drugs available for treatment (referred to as actionable). The selected 61 target gene exons are those that have been reported to be the junctional driver gene exon in such fusion genes in the literature and in clinical databases. 152 primers (Table 3) were designed to cover these targets and are used in a multiplexed primer pool to amplify target exons from ds cDNA with anchored linear PCR using the Qiagen Multiplex PCR kit.

[0147] Table 3: Sequences of Outer Primers presented in 5’ to 3’ orientation.

[0148] The multiplex primer mix covering these exons was made at 1 uM and a final concentration of 0.1 uM was used in a 50ul reaction with 11 ,25ul of template cDNA, and 30ul of the Qiagen Multiplex 2x mastermix. Thermalcycler program: 95C 15 minutes, 8 cycles of (94C 30 seconds, 60C 3 minutes, 72C 30 seconds), 4C hold. Enrichment product was purified with a 1 ,8X Ampure XP SPRI purification and elution in 20uL EB.

[0149] Adapter Ligation: Resulting ssDNA amplicons were used as ligation phosphoacceptor template molecules in a ligation with a ssDNA phosphodonor adapter coordinated by a splint RNA in a DNA:RNA hybrid duplex (Figure 3) to permit PCR-based amplification of the linear amplicons using NEB SplintR ligase. The ligation process involved three steps and oligos listed in Table 4.

[0150] Table 4: Sequences of Ligation Oligos presented in 5’ to 3’ orientation. 75Phos / ’ refers to a 5’ phosphate modification. 73SpC3 / ’ refers to a 3’ C3 spacer modification, ‘r’ indicates subsequent base is a ribonucleotide while ‘N’ indicates random incorporation of any of the four standard bases at the given position.

[0151] Reaction 1: Annealing of ssDNA phosphodonor adapter and RNA splint (Figure 3): 95u I of 10OuM adapter and 95 u I of 10OuM splint were combined with 10ul of annealing buffer (10mM Tris pH 8, 1 mM EDTA, 1 M NaCI). Thermalcycler program: 85C 2 minutes, 80 cycles of 30 seconds with a reduction of -1C per cycle, 4C hold.

[0152] Reaction 2: Annealing of ssDNA phosphoacceptor amplicon with antiprimers pool (Table 5) to prevent adapter ligation to any unextended outer primer not removed by SPRI clean-up: 19ul of enriched amplicon was combined with, 0.125uL of 2uM antiprimers (0.25pmol), 1.15uL annealing buffer, and 2.725ul of water.Thermalcycler program: 95C for 2 minutes, 70 cycles of 30 seconds with a reduction of -1C per cycle, 25C hold.

[0153] Table 5: Sequences of Antiprimers presented in 5’ to 3’ orientation.

[0154] Combination of reactions 1 and 2: The entire reaction 2 product (23uL) is combined with 7uL of a ligation mastermix (1 uL of 50uM annealed splint:adapter, 3uL of 10x NEB splintR reaction buffer, 2.25uL EB, 0.75uL NEB splintR ligase). Thermalcycler program: 25C 8 hours, 9 cycles of (6 minutes, -1C per cycle), 16C 30 minutes, 65C 2 minutes, 4C hold.

[0155] Ligation product was purified with a 1 ,8X Ampure XP SPRI purification and elution in 20uL EB.

[0156] Nested Target PCR Amplification: The ligation product was amplified using a general amplification primer i7R with sequence complementarity to the ligated adaptor molecule (Nextera i7_R: GTCTCGTGGGCTCGGAG, SEQ IDNO: 3), and a nested gene-specific primer pool which included 150 gene specific primers (GSP) (Table 6) that are internal to the first gene specific linear amplification primer using KAPA HiFi hotstart readymix.

[0157] Table 6: Sequences of Inner Primers presented in 5’ to 3’ orientation. Primers consist of a common 5’ Nextera XT universal amplification sequence followed by a target exon specific GSP sequence.

[0158] The nested primer mix was made at 2uM, and a final concentration of 0.1 uM nested primer mix and 0.1 uM Nextera i7 R primer was used in a 50ul reaction with 19uL ligation product. Thermalcycler program: 95C 15 minutes, 8 cycles of (94C 30 seconds, 60C 3 minutes, 72C 30 seconds), 4C hold. Enrichment product is purified with a 1 ,8X Ampure XP SPRI purification and elution in 20uL EB.

[0159] Sequencing Library Construction: A sequencing library is generated by the PCR-based addition of indexing and clustering NGS tails using the Nextera XT system (Illumina) and IDT for Illumina UDI adapters (Illumina) according to manufacturer instructions. 7.5ul amplification product is combined with 12.5ul KAPA HiFi PCR mix, and 5ul adapter mix in a 25ul reaction. Thermalcycler program: 95C 30 seconds, 15 cycles (95C 30 seconds, 55C 30 seconds, 72C 30 seconds), hold 4C. Library is purified with a 0.7X Ampure XP SPRI purification and elution in 30uL EB.

[0160] Sequencing Library Normalization, Dilution, and Sequencing: The libraries were quantitated with a qubit BR DNA kit (Thermofisher) according to manufacturer instructions and individually normalized to 20nM. 5uL of each library isthen pooled, and the pool was diluted to 4nM before denaturing and diluting according to Illumina MiSeq protocol. The final library was sequenced at 12.5pM with 5% PhiX using an Illumina MiSeq instrument with a 300 cycle cartridge.

[0161] As a demonstration, we performed ligation in the absence of individual components to demonstrate molecular shifts in size of the apparent 400bp ssDNA template generated from a pool of 78 synthetic dsDNA “geneblocks” (Table 7) with an average length of 500bp and with each construct encoding a fusion gene construct involving one of the targeted gene exons (Figure 12).

[0162] Table 7: Sequences of Geneblocks presented in 5’ to 3’ orientation.

[0163] As indicated in Figure 12A, antiprimers (reverse complement oligos of the GSP molecules) were included to convert excess unextended gene specific primer into dsDNA duplexes not available for ligation. Figures 12B-G show visualization of 1 ng of components with Agilent Bioanalyzer High Sensitivity DNA Chip, as follows (B) Linear ssDNA template (C) With ligation and antiprimer (D) With ligation without antiprimer (E) Without ligation (no adaptersplint) (F) Without ligation (no ligase) (G) dsDNA template. electrophoresis as electropherogram and gel image composite. Ligation was expected for Figures 12C and 12D and not for 12E and 12F.

[0164] Because the bioanalyzer system uses a dsDNA ladder and ligation products are ssDNA, partial duplexes, hybrid duplexes and dsDNA duplexes, migration sizing was not accurate. We observed the initial ssDNA template and two shifted populations representing ligation negative and ligation positive reactions. In the absence of splintadapter hybrid duplex or ligase (ligation negative) we observed one shift. In the presence of all components and in the absence of antiprimers, which are annealed to the linear amplicon to prevent ligation of unextended primer, we observed a second shift (ligation positive). This is as expected as both the splintadapter and ligase are necessary for ligation while the antiprimers are not and therefore their absence does not affect ligation at the phosphoacceptor 3’ end however because they form partial duplexes with the phosphoacceptor molecules at their 5’ ends they do affect migration.

[0165] Following ligation of phosphoacceptor linear ssDNA amplicon product with ligase to phosphodonor adapter via random splinting, next generation sequencing libraries were generated from ligation product by performing PCR with a primer specific to the phosphodonor adapter (Nextera i7_R: 5’-GTCTCGTGGGCTCGGAG-3’, SEQ ID NO: 3) and with nested primers internal to the first GSP primer set used to generate the ssDNA ligation template (Table 6) (Figure 13). These primers have Illumina Nextera platform tails that were subsequently used to add indexes and add Illumina p5 and p7 platforms.

[0166] Following ligation of geneblock derived linear amplicon, next generation sequencing libraries were generated from ligation product (Figure 14). When libraries were generated from reactions lacking initial dsDNA template, we observed no difference in library mass yield when ligation components were all present or individually omitted. When geneblock template was used, we observed a 15-35x difference between ligation negative and ligation positive reactions. In ligation negative reactions non-specific amplification was observed while in ligation positive reactions specific product was observed at 450bp.

[0167] To demonstrate the applicability of this method to an extensively diverse population of endogenous cellular RNA, an Illumina NGS library was generated from universal normal human RNA (Figure 15). Total RNA was randomly converted into double stranded cDNA which was then subjected to linear amplification by PCR to generate linear ssDNA template for ligation. When the ligation product was amplified by PCR, non-specific product < 300bp was observed for a library from ligation product made without template as well as from product in the absence of ligase. Specific library of size 300-1 OOObp was generated in the presence of ligase indicative of a diverse population of amplicon with variable lengths. Library yield was 3fold that of the ligation negative and template negative libraries.

[0168] When ligation product was converted into Illumina NGS library and sequenced via Illumina paired end 300 cycle miseq sequencing with ~1 M reads per sample, the diversity of individual template molecules could be assessed by the counting of unique amplicon end positions represented in the library since template RNA molecules exist as a randomly fragmented population (Figure 16). When sequenced, the end position is the start of the illumina P7 read while the GSP sequence is the start of the illumina P5 read. To assess diversity, we examined TUBB gene amplicons observed in the library generated from universal normal human RNA. The TUBB cDNA transcript was amplified with a single gene specific primer 608bp from the start of the TUBB transcript meaning that amplicons of up to 608nt could be theoretically observed. Due to degradation of the 5’ end due to RNA instability and the inability of random cDNA priming to capture the entire 5’ end of thetranscript observation of the entire transcript is however not expected. Consequently, we observed end sites of up to 529nt from the start of the GSP and saw that longer transcripts are more relatively more rare than short transcripts as predicted due to RNA instability.

[0169] To quantify the diversity of initial template molecules incorporated into the library, NGS libraries were prepared from multiple templates: no template control (NTC), amsbio universal normal human RNA, seraseq v4 RNA fusion mix, and a blend of amsbio universal normal human RNA and synthetic geneblocks (Table 8).

[0170] Table 8: Read count and template diversity from ligation product libraries

[0171] Libraries were sequenced via Illumina paired end 300 cycle miseq sequencing with ~1 M reads per sample. Gene specific primers were used in the linear amplification process to amplify regions from 150 exons within 21 genes. As expected, no significant number of on-target reads were observed when library was generated without template. In contrast, 244,954-533,337 on-target reads were generated from template positive samples. This corresponded to 6224-15,086 unique on-target ends suggesting at least this many unique template molecules were incorporated into the library. These reads were associated with 20-21 of the 21 expected genes and 103-148 of the expected 150 exons. Coverage is not expected for all targets in all samples due to differences in gene expression and input amounts.

[0172] Accordingly, application of the fusion detection methods permitted identification of structural rearrangements in reads from the previously described sample libraries after Illumina paired end 300 cycle miseq sequencing with ~1 M reads per sample (Table 9).

[0173] Table 9: Read count, template diversity, fusion calls, fusion read depths, and fusion unique ends from ligation product libraries

[0174] The NTC sample is template negative and fusion negative and consequently no significant number of reads and no fusions were detected. The amsbio universal normal human RNA sample is template positive with no fusions expected and this was observed. The seraseq v4 RNA fusion mix is template positive and contains 18 expected fusions which were all detected. The amsbio universal normal human RNA sample spiked with 10pg of 78 fusion geneblocks resulted in detection of all detectable 56 fusions based on the target GSP sequences used.

[0175] The resulting sequencing library, with inserts derived from the set of PCR amplified adapter- ligated DNA molecules, is then subjected to computational analysis that includes filtering and error correction steps for the purposes of identifying among this set, any target gene with atypical exon configuration, or target gene which has been fused with any other DNA sequence such as a partner gene in a fusion configuration. Such genomic rearrangements or fusion events are clinically actionable with targeted inhibitors. In the case of cancer for example, targeted inhibitors can be administered to cancer patients when such fusion events are identified in patient samples. This combinational analysis is now described.

[0176] The above-mentioned computational analysis is referred to as a fusion pipeline which is programmed and executed locally on one or more computer systems, or in the cloud.

[0177] Input Files

[0178] Input files stored on some medium are provided to the computing system. These input files can include data stored or presented in any conventionally known format. Alternatively, if conventional formats are insufficient for the fusion pipeline, then custom formats can be developed with the fusion pipeline algorithmconfigured to recognize these custom formats. The required files for the fusion pipeline of the present example include:

[0179] - A reference genome FASTA file, assembly GRCh38.p13, openly available at the Ensembl FTP server https: / / uswest.ensembLorg / Homo sapiens / lnfo / lndex.

[0180] - An annotations file in TSV format, with all the CDS exons and corresponding transcript annotations (Ensembl 99, obtainable at Ensembl Biomart -

[0181] The assay manifest file, with information about all gene targets and corresponding gene-specific primers. This can be the information previously described in the target enrichment section above.

[0182] - The test samples FASTQ files, the output of the Illumina MiSeq sequencer as described in the Laboratory Methods section, containing the cDNA sequences where fusion events will be called. This information is obtained after using the adapter / amplicon / primer PCR methods that were previously described. For each sample, the sequencer outputs two files, one with the R1 reads, corresponding to reads starting at the gene-specific primer, and another with the R2 reads, starting from the adapter end. R1 and R2 reads are referred to as read 1 and read 2 on Figure 2, in the Laboratory Methods section. Each R1 read has a corresponding R2 read and are called a read pair.

[0183] Read Filtering

[0184] The first step in the pipeline is to process reads from the FASTQ files, filtering out-of-target / misprimed reads and also generating QC metrics.

[0185] For each R1 / R2 read pair, since it is known that the R1 read starts at a gene-specific primer on a driver gene (the known sequence, Figure 2), the first step is to identify which primer is present at the start of the R1 read, with maximum hamming distance of 1 from all the primer sequences in the assay manifest file. If there is a primer satisfying this distance condition, the read pair is labelled with the corresponding primer, otherwise the read pair is removed.

[0186] Then, for all labelled read pairs, the start of the R1 read is aligned with the expected target sequence (primer + 20 bps, up to exon boundary, again theknown sequence of Figure 2). The alignment distance dj between a read / and the expected target sequence j is defined aswhere I is the length and G,j is the global alignment score of the compared sequences, with score 2 for match, -1 for mismatch and no penalty for gaps, obtained with a dynamic programming alignment using the open-source Bio.pairwise2 module from Biopython. If the primer binds to the known sequence region, as expected, this distance should usually be zero, so larger values point to an off-target sequence (when the primer binds to a different region) and should be filtered out. But sequencing errors can also cause the distance to increase, even when the primer correctly amplified the known region, so the challenge is to select a threshold that is small enough such that most off-target sequences are filtered out, without removing good sequences that have sequencing error. After several tests, the distance threshold with best results was found to be 1 . Therefore, all read pairs where> 1 are filtered out, and the remaining are saved in an intermediate file to be used as an input in the fusion calling step.

[0187] Fusion Calling

[0188] The fusion calling step takes as an input the filtered reads file generated in the read filtering step above, in addition to the reference genome and annotations file. Each of the phases of this step are described in the sections below.

[0189] K-mer Indexing

[0190] In order to identify the origin of a particular cDNA sequence, a mapping step is performed, where a given read sequence from a test sample is mapped to a reference sequence (in this case, the human genome), so the step of finding genome mutations can be done. This is a common step in any variant calling pipeline, and different methods have been proposed. In the described pipeline according to the present embodiments, a k-mer mapping algorithm is used, and a few definitions will be needed to describe this mapping step. It should be noted that this fusion pipeline and any of its described constituent steps are not limited to human genome analysis. In alternate embodiments, the fusion pipeline can be configured to analyze any animal genome, provided a reference sequence is available.

[0191] A k-mer is a substring of length k of a given DNA sequence. For a given k-mer, the canonical representation, or simply canonical k-mer, is the lexicographically smaller of the k-mer itself and its reverse complement. For instance, for the 5-mer TAGCT, the reverse complement is AGCTA and therefore AGCTA is the canonical k-mer. This means that a k-mer and its reverse complement have the same canonical representation.

[0192] A k-mer index is a structure where information about a given k-mer from a (usually very large) set of k-mers can be accessed. This information is commonly the list of sequences and locations where this k-mer is present, such a reference genome or transcriptome.

[0193] In this pipeline, the k-mer index I is a simple hash table where the keys are all canonical k-mers found in the reference genome, and the information stored is the set of tuples, with exons that have the k-mer and the coordinate of the k- mer hit within the exon. The default k-mer length is 17, however the entire process can be repeated with a different value of k. The default value of 17 was determined based on tests with different k-mer sizes. In similar projects, either a default is selected based on performance tests, or a range of different values can be used, and all results are compared. In one test case, a k value of 17 was optimal, but k of 19 or 21 are also usable with similar results. I[k] is defined as the k-mer hits of a k-mer k: l[k] = {(exon_1 , pos_1), (exon_2,pos_2), ...}

[0194] For instance, in the example from Figure 5 (a), l[GCAGGTT] = {(ENSE001 , 1)} and l[TGGCAGG] = {(ENSE001 , 12)}. ENSE001 is an example exon identifier, based on the Ensembl nomenclature, where all the exon identifiers have the “ENSE” prefix.

[0195] Alternative structures have been proposed in the literature to store k- mer indexes in a more efficient way, especially with respect to memory usage reduction and might be used on future versions of the k-mer index.

[0196] The k-mer index is created on step 2 on Figure 4. It is kept in memory during the execution, but other approaches could be used in order to minimize the memory usage, such as using a fast database suchthat would offer a trade-off between speed and memory.

[0197] After its creation, the k-mer index is stored on a network or some removable storage medium in a serialized form using picklewhich means that it does not need to be created for each new run of the pipeline, if there are no changes in the reference genome or annotation file, and can instead be directly loaded from the pickle file.

[0198] Read Pseudo-alignment

[0199] This step (Figure 4, step 3) uses the k-mer index to map reads to reference genome exons. The reference genome is the linear DNA sequence with all chromosomes. Then with the k-mer index the system can find the genome location of a read by analyzing the k-mer patterns of the read. The aim is to identify the corresponding exon regions of the sequencing read (Figure 5 (b)), to then find reads that map to two different genes (the driver gene in the known sequence, and a different partner gene in the unknown sequence, Figure 2), indicating a possible fusion supporting read. The following subsections describe each step in the pseudoalignment phase.

[0200] K-mer mapping and finding exon candidates

[0201] Given a k-mer from a sequencing read, a k-mer hit is the position of this k-mer in the read and the mapping information of this k-mer from the k-mer index. If a k-mer is not present in the index, the corresponding mapping information is called a miss.

[0202] For instance, the k-mer hit of a k-mer k_1 on position 1 of a read is a tuple (1, l[k_1]), where l[k_1] is the k-mer index information of this k-mer (or miss if the k-mer is not present in the index). Therefore, repeating this process for all k-mers in a read results in a vector of k-mer hits, representing the original read:(1 , l[k_1]), (2, l[k_2]) (n, l[k_n])Where n is the number of k-mers on the read, and k_1, k_n are all the k-mer in the read.

[0203] Consecutive hits for the same unique exon match are then aggregated into a single element called exon candidate ej, which is defined as a tuple of three elements: the exon, the range of the hit on the read, and the range of the mapped region on the exon: (e_1 , read_st— >read_end, map_st— >map_end). Consecutive misses are also merged, but instead of a mapped range, the length is stored instead.The combined misses are called a missed block. Exon candidates with small range are filtered, that is, if read_end - read_st < min_threshold. The default min_threshold is 24 - k + 1 , resulting in a minimum mapped requirement of 24 bps in the present embodiment.

[0204] At the end of the process, a read is then represented as a mapping vector r_i composed of consecutive exon candidates and misses.

[0205] For instance, r_i = [(e_1 , 1 -^41 , 50— >90), (miss, 42— >61 , 20), (e_2, 62— >85, 1 — >24)] means that from the base pairs 1 to 41 , the read maps to region 50- 90 on exon e_1 , then the next 20 base pairs on the range 42-61 do not map to anything on the k-mer index, and finally the base pairs from 62 to 85 map to region 1 to 24 on exon e_2.

[0206] Sequencing error / SNP correction

[0207] This error correcting step is optional if sequencing errors are at a minimum or nonexistent. This error correcting step is applied to any read vector from the previous step that has a missed block, such that the corrected version is retained in the pool of exon candidates. Read errors or SNPs might cause the k-mer mapping step to miss some hits. When this happens, it is possible to correct the resulting mapping if there are exon candidates “around” the missed k-mers, if the mapped coordinates agree. Specifically, consider a missed block with exon candidate block with the same exon before and after it:(e_1 , r_1 — >r_2, m_1 — >m_2), (miss, len), (e_1 , r_3— >r_4, m_3— >m_4)If r_2 + len + 1 = r_3 (or similarly, m_2 + len + 1 = m_3), then we can assume that the missed block was caused by an error since the surrounding areas agree.

[0208] For instance, in the example shown on Figure 6, the initial mapping results in rj = [(e_1 , 1->3, 51->53), (miss, 4— »10, 7), (e_1 ,11 ^12, 61 ^62)]. Since 53 + 7 + 1 = 61 , the miss can be removed and the two exon candidates are merged, resulting in rj = (e_1 , 1— >12, 51 — >62)].

[0209] Resolving mapping ambiguities

[0210] This step is applied against all the valid exon candidates resulting from the previous steps above. After the k-mer mapping step, some mapped exon candidates might be overlapping, either from mismapped regions due to sequencesimilarity between different regions, or different exon annotations for the same genomic region, among other reasons. Those ambiguities are resolved, resulting in a mapping without overlapping exon candidates. Following is an example of how a mapped ambiguity is corrected.

[0211] For instance, if rj = [(e_1 , 1 ^21 , 51^71), (e_2, 7^15, 62->70)], meaning that region 1-21 of this read mapped to e_1 from 51-71 , but also region 7- 15 mapped to a different exon e_2, from 62-70. Therefore, there is an overlap (in this case the mapping of e_2 is completely contained in e_1 , but not always). The e_2 mapping is removed, since it is smaller than e_1.

[0212] For each overlapping pair of exon candidates, the following checks are made until one of the conditions is met:1. If the exon is the same on both candidates, they are merged, adjusting the mapping coordinates accordingly.2. If one of the candidates has a much smaller mapping range than the other, it is removed. This is because a larger mapped area indicates a higher confidence in the accuracy of the mapping. The chance of a spurious match decreases rapidly with increasing mapped length.3. If one of the candidates has an exon from an assay targeted gene, the other candidate is removed. It is assumed that sequences come from the targeted genes, so if there is ambiguity in a match, because the sequence is similar to different parts of the genome, this is resolved in favour of the targeted gene match, since it is the expected sequence to see.4. If one of the candidates matches a list of preferred candidates, the other candidate is removed. This is useful when mapping the R2 read, using the exon candidates from R1 as the preferred candidate list.5. If both candidate exons are from the same gene, choose the candidate with highest priority, calculated from sorting the corresponding transcripts for each exon in descending order of length.6. If all previous tests do not solve the ambiguity, delete the candidate with the largest exon gap, and in case of a tie, keep the candidate with lexicographically smaller ID, to be consistent.

[0213] Amplicon End Estimation

[0214] Since the assay is based on an anchored linear amplification step, different sequenced amplicons will have different lengths, and therefore different endpoints (or starting R2 sequences). This step is a data extraction of information of the various amplicons to determine their uniqueness based on its endpoint and length. The genomic coordinate of the amplicon end of the current read pair being processed can be estimated by looking at the R2 read mapping. In this mapping, the end coordinate estimate is obtained by adding the last exon candidate mapping coordinate to the length of a missed region (if present, to account for errors at the end of the read). On the other hand, if there is no missed region then nothing is added. Using the previous example, if r_i = [(e_1 , 1 — >41 , 50— >90), (miss, 42— >61 , 20), (e_2, 62— >85, 1 — >24)], there is no missed region at the end, so nothing is added, and 24 is the end coordinate. But suppose that rj = [(e_1 , 1 -^41 , 50— >90), (miss, 42— >61 , 20), (e_2, 62— >85, 1 — >24), (miss, 86 -> 90, 5)]. Then 5 is added to 24, giving 29 as the end coordinate.

[0215] This will be important to estimate the number of unique amplicons on the list of supporting reads for all fusions. It is noted that all amplicons have the same starting point, since they all start at primer. However, the endpoint is different due to the application process.

[0216] Mate-pair Merging

[0217] After the k-mer mapping steps described above have been performed on both reads of a read pair, a merging algorithm tries to combine both reads into a single mapping vector. The algorithm merges all exon candidates from both mapping into a single vector (excluding mismatches) and sort it by mapping coordinates ( After that, apply the algorithm for solving mapping ambiguities presented in the section above. This will merge concordant exons on both mapping and fix possible discrepancies. The result will be a single mapping vector, per read pair, that will then be used to find fusion events.

[0218] Figure 7 graphically shows an example of this process, for a read pair mapping to three different exons, e_1 , e_2 and e_3. The resulting mapping for both reads are r_1 = (e_1 , 1->20, 21 — >40), (miss, 21 ^40, 20), (e_2, 41 ^80, 1->40), (miss, 81 ->100, 20), (e_3, 101 ^120, 1 ->20) and r_2 = (e_2, 1->30, 11 ->40), (miss, 31 ^50, 20), (e_3, 51 ^90, 1->40), (miss, 91 — >120, 30).After combining and sorting by genomic coordinates, the result is the mappingThen, the merging process is applied, resulting in r_m = (e_1 , 1 — >20, 21— >40), (e_2, 1->80, 1— >40), (e_3, 51— >120, 1->40).

[0219] Junction Candidates

[0220] After the mate-pair merging, the next step is to find possible fusion junctions on the resulting mapping vector.

[0221] A junction candidate is defined as a pair of consecutive exon candidates on a given read mapping, with a miss region in between. If each exon is from a different gene (a known driver gene and an unknown partner gene), it will be called a fusion junction candidate, otherwise it is a normal junction candidate. To find junction candidates, the algorithm scans the mapping for consecutive exon candidates, with some additional conditions described below.

[0222] Consider that for a “perfect” fusion, where the entire sequence of the driver gene exon is immediately followed by the entire sequence of a exon from a different partner gene, a particular “fusion signature” in a mapping vector is expected as shown in Figure 8. In Figure 8, the 5-mers of exon 1 (shown below the box labelled Exon 1 , representing the driver gene) are followed by four 5-mer misses (shown above the boxes labelled Exon 1 and Exon 2), then by the 5-mers of Exon 2 (shown below the box labelled Exon 2, representing the partner gene).

[0223] where the length of the missed block is exactly k-1, and k is the k-mer length. This happens because the k-mers that cover the junction of the two exons are not naturally present in the DNA and therefore are not also in the k-mer index (except for spurious hits that are most likely filtered in the previous steps). There is exactly k- 1 k-mers covering the junction, which results in the missed block with length k-1 in the mapping vector.

[0224] In addition, if the exons are completely present, the junction point happens exactly between the extremities closest to each other of both exons. This is referred to as a perfect signature as shown in Figure 9a. In this example, m_2 is either exon e_1 length or 1 (the latter meaning that the mapping is in reverseorientation), and similarly m_3 should also be either exon e_2 length or 1 . The exon gaps g_1 and g_2 is defined as the distance between m_2 and m_3 to their respective exon extremities. Also, the k-mer gap g_k is defined as the difference between the length of the miss block and (k-1), in the junction signature. The junction signature gap (g) is then defined as the tuple (g_1 , g_k, g_2). In a “perfect” signature, g = (0,0,0).

[0225] Junction signature gap in the presence of sequencing errors

[0226] If sequencing errors or SNVs happen close to the junction point, the previously discussed error correction techniques will not work, since there are no exon candidates for the same exon on both sides of the error. More specifically, because the error correction is based on having good hits around the error, on both sides. If the error is too close to the junction (less than the k-mer size), then there will be no good hits on one of the sides, closest to the junction point. Also, the k-mer misses caused by the error will result in changes in the fusion signature gap, as can be seen on Figure 9a-9d. Figure 9a is an example showing no errors, or a perfect signature, with g = (0,0,0). Figure 9b and 9c shows examples of one error 2 bps from the junction point, causing two k-mers to miss, increasing both the exon gap g_1 and the k-mer gap by two. Figure 9d is an example showing e.rrors on both exons, causing g_1 = 5, g_2 = 2, and the k-mer gap g_k = g_1 + g_2 = 7.

[0227] In the presence of errors, the exon gaps g_1 and g_2 increase by the same amount as the k-mer gap g_k, leading to the identity g_k = g_1 + g_2. This fact points to a gap correction algorithm, where g_k can be reduced by one, also reducing either of the exon gaps (the largest), until either g_k is zero, or both g_1 and g_2 is zero. The gap correction algorithm is expressed as follows below, where sub-parts a and b are iterative until the listed conditions are met.Gap correction algorithm:1. While g_k > 0 a. If g_1 > 0 then i. g_k = g_k - 1 ii. g_1 = g_1 - 1 b. Else If g_2 > 0 then i. g_k = g_k -1 ii. g_2 = g_2 - 1

[0228] Clustering fusion junctions and recruiting supporting reads

[0229] For all junction candidates found in a read pair, if there is a single fusion junction candidate where the k-mer gap is smaller than a given threshold (default is k + 1), then the read pair will be added to a list of supporting reads for this fusion junction.

[0230] It should be noted that the previously described gap correction reduces the gap in some cases where there are errors close to the junction, but it does not remove all gaps. Gaps can still be large if there are different sections of the read mapping to different parts of the genomes, and such gaps cannot be corrected, but if they are still large it mostly means a mapping error or random mapping and will not be considered a fusion. This mostly happens due to sequence similarities / repeats on different parts of the genome, meaning that a read can map to different regions but does not necessarily means that the read comes from that region, is just a mismatch due to seq similarity.

[0231] If a read has no fusion junction candidates but it has normal junction candidates, it will be counted as a normal junction forthat particular gene. This is useful to estimate the ratio between the number of fusion supporting reads and normal reads for all fusion calls.

[0232] After all read pairs are processed, all fusion junction candidates are clustered based on the genomic coordinates of the junction point, with a tolerance of 50 bps, meaning that junctions that are less than 50 bps apart are considered to be the same fusion.

[0233] Fusion Candidate Calling

[0234] After the junction clustering step, a list of fusion candidates for the sample results, each with a number of supporting read pairs. For each fusion candidate, several metrics are calculated based on the supporting read pairs, such as the read depth, the number of unique amplicon ends, the estimated length of the fusion transcript, the ratio between fusion and normal reads, among other metrics. The estimated length is based on the average length of the read mapping length of all supporting reads for the fusion.

[0235] Local assembly of fusion transcripts

[0236] For each called fusion, all supporting reads for that fusion are assembled using de-Bruijn graph-based assembly. After the assembly, all reads of the sample are mapped against all assembled fusion transcripts. This is done as follows:1. For each fusion: a. a de-Bruijn graph (https: / / en.wikipedia.org / wiki / De Brush graph) with the k-mers of all supporting reads is built. b. The graph is simplified with bubble correction and tip removal algorithms, to account for sequencing errors. c. The graph is traversed, starting with the node that has the k-mer of the gene-specific primer that captured the fusion, and a number of contigs is extracted. d. The fusion calling algorithm is rerun on the output contigs. The longest contig that has the expected fusion junction is considered the correct assembled fusion transcript.2. A reference file with all fusion transcripts is created, and the original FASTQ files from the sample are mapped to this reference file using BWA (http: / / bio- bwa rce io rg e , net) .3. New metrics are added to the fusion calls, such as the assembled transcript length, the number of total mapped reads and unique mapped reads, and the estimated number of unique transcript ends from the mapping output.

[0237] Confidence Score Assignment

[0238] The fusion calling algorithm is lenient since it looks for reads with mapping regions in two different genes. Although there are some filtering steps, there is at this point a high number of artifact calls. Therefore, all called fusions will be scored based on the corresponding metrics and depending on the resulting score will be either filtered and discarded, or receive a confidence label, from Low, Medium or High. An overview of the confidence scoring method can be seen on Figure 10.1 . Artifact filter: The first step checks if the orientation of both fusion partners is concordant, and if exactly one of the partners is a gene targeted by the assay. If either of the conditions is not satisfied, the fusion is filtered as an artifact. Fusion junction candidates in the merged mapping vectors having a pair of consecutive exon blocks from different genes are called fusion partner genes.2. Supporting fusions: Even in a normal cell, because of alternative splicing, the same gene can have different transcripts (isoforms), meaning different variants of the same gene that skip / include different exons. This means that the same fusion, caused by a breakpoint in the genome, can also have different variants with the same two partner genes, but between different exons. For example, if it is seen in the same patient a fusion EML4;7 - ALK;20 (between EML4 at exon 7 and ALK exon 20) and a fusion EML4;6 - ALK; 20. Instead of considering two different fusions, the one with the most supporting reads is considered the main fusion, and all others are supporting fusions, meaning that the supporting fusion reads are added to the main fusion. This in turn can improve the main fusion score.3. Metrics filter: Several metrics are checked each with a minimum threshold value that the fusion has to pass. Failing on any of the metrics will label the call as filtered.4. High confidence: If a fusion call passed the previous steps and has read depth and number of unique reads higher that a given threshold, it will be labelled high confidence.5. The remaining fusions will be processed by a scoring function that will be explained below. Depending on the score, a fusion can either be considered medium or low confidence, or be filtered.

[0239] Fusion scoring function

[0240] For a given fusion j, and a set of metrics, the fusion score s(j) is given by metricswhere m_ij is the value of the metric / for the fusion j, wj is the weight associated with metric / , and p_\ is the expected mean of metric / , calculated from a control dataset with known fusions.

[0241] The set of metrics can be chosen from several different metrics, with varying results. Best results in our tests were obtained by using four metrics: read depth, number of unique reads, fusion transcript length, and the fusion / normal readdepth ratio. The fusion / normal read depth ratio is the total number of normal junctions for the same gene as the target gene for the fusion.

[0242] By looking at the scores of fusions in a training set, two thresholds t_ 1 and t_2 are determined, such that fusions with score > t_1 are labelled medium confidence, while fusions with score < t_2 are filtered and discarded. The remaining fusions receive a low confidence label. Figure 11 shows the scores for fusion calls made on a control dataset with known fusions. The thresholds were defined as t_1 = 0 and t_2 = -1. Scores that are categorized as medium and low confidence are subject to review for further consideration of actions to be taken, including additional analysis. For each confidence level, different performance metrics can be determined, in terms of sensitivity and specificity. For example, low confidence means that there is a significant chance that the fusion might be an artifact, so an orthogonal validation method should be used for confirmation. Conversely, high confidence calls are >99% guaranteed to be correct, so no validation should be needed. Medium is an in-between confidence, so it’s up to the discretion of who is analyzing the results to determine what action to take. If the result is expected, such as a common fusion for a specific cancer type, then no confirmation might be needed.

[0243] Following is a summary of the fusion pipeline method. In a first aspect, the fusion pipeline method is a computer implemented method for identifying from a plurality of read pairs R1 and R2 of a sequenced sample with inserts derived from a set of PCR amplified adapter ligated DNA molecules, any target gene fused with any other DNA sequence which are potentially clinically actionable. The method includes generating a k-mer index of predetermined length k from a reference genome or a transcriptome; mapping each sequencing read from the plurality of read pairs to exons of the reference genome with the k-mer index to identify corresponding exon regions of the sequencing read, and generating for each sequencing read a mapping vector including one or more blocks, where each block is either an exon block defined by a range of consecutive k-mer matches and their corresponding locations on an exon of the reference genome, or a mismatch block defined by a length of consecutive k-mer mis-matches;identifying exon blocks with overlapping consecutive k-mer matches in each of the mapping vectors, and merging the same exon blocks with adjustment of mapping coordinates or retaining the exon block having predetermined high priority criteria to resolve mapping ambiguities; combining both the mapping vectors corresponding to each read pair and sorting by genomic coordinates to provide intermediate mapping vectors, and merging any exon blocks with the same respective exons with adjustment of mapping coordinates to generate merged mapping vectors; identifying fusion junction candidates in the merged mapping vectors by searching for a pair of consecutive exon blocks from different genes, called fusion partner genes, with a mismatch block in between, and recording the corresponding read pairs for each of the identified fusion junction candidates; clustering all the fusion junction candidates based on genomic coordinates of junction points of the pairs of consecutive exons from different genes where junctions closer than a predetermined number of base pairs are considered to be the same fusion junction candidate, and calculating a first predetermined set of metrics based on the read pairs including read depth; assembling a fusion transcript for each fusion junction candidate cluster using the corresponding recorded read pairs to generate assembled transcripts, and calculating a second predetermined set of metrics from the assembled transcripts; mapping the plurality of read pairs to the fusion transcripts, and calculating a third predetermined set of metrics based on the mapping output including a number of mapped reads; and categorizing any fusion transcript as being a high confidence reportable fusion when the orientation of both fusion partners is concordant, exactly one of the partners is a gene targeted by the assay, minimum thresholds for specific metrics of the first, second and third predetermined set of metrics are met, andthe number of mapped reads and the read depth exceed predetermined thresholds.

[0244] In a first embodiment of the first aspect, exon blocks having less than a minimum mapped bps are filtered out, and the minimum mapped bps is 24 bps.

[0245] In a second embodiment of the first aspect, each of the mapping vectors is represented by the expression (e_n, read_st— >read_end, map_st— >map_end) for each exon block and (miss, miss_st— >miss_end, len) for each mismatch block, where e_n is the exon unique identifier, read_st and miss_st are the starting base pair on the sequencing read, read_end and miss_end are the ending base pair on the sequencing read, map_st is the starting mapped region on the exon for read_st, map_end is the ending mapped region on the exon for read_end, miss is a status identifier for a mis-match, and len is the bps length of the mismatch block.

[0246] In a third embodiment of the first aspect, the method can include error correcting the mapping vectors to remove mismatch blocks interposed between two exon blocks with the same exon. In this embodiment, error correcting includes removing the length of consecutive k-mer mis-matches of the mismatch block interposed between a first range of consecutive k-mer matches of a first exon block and a second range of consecutive k-mer matches of a second exon block of the same exon, when a continuation of the first range by the length of consecutive k-mer is sequentially followed by a beginning of the second range. More specifically, error correcting to remove a mismatch block between two consecutive exon blocks with the same exon includes determining when read_end + len + 1 for the first exon block = read_st for the subsequent exon block to remove the mismatch block. Alternately, error correcting to remove the length of consecutive k-mer mis-matches includes determining when map_end + len + 1 for a first exon = map_st for a subsequent exon to remove the length of consecutive k-mer mis-matches.

[0247] In the second embodiment, identifying fusion junction candidates includes executing an error correction applied on k-mer gaps g_1 , g_2 and g_k of the fusion junction candidates when g_k > 0, where g_1 is a distance between map_end of the first exon block to the corresponding exon extremity,g_2 is a distance between map_st of the second exon block to the corresponding exon extremity, g_k is the difference between the k-mer length k and the mismatch block len, and the k-mer gap correction steps include while g_k > 0 a. if g_1 > O then i. g_k = g_k - 1H. g_i = g - 1 b. else If g_2 > 0 then i. g_k = g_k -1 ii. g_2 = g_2 - 1

[0248] Alternately for the second embodiment, the predetermined high priority criteria include an exon block with a larger mapping range than the other, an exon block where its exon comes from an assay targeted gene, an exon block matching a list of preferred candidates, or when both exons are from the same gene, selecting the exon with highest priority, calculated from sorting the corresponding reference genome transcripts for each exon in descending order of length.

[0249] If none of the predetermined high priority criteria is met, then the exon block having the largest exon gap is removed, where an exon gap is a distance between map_end of a first exon and map_st of a subsequent exon to their respective exon extremities.

[0250] In this embodiment, assembling the fusion transcript includes the following steps for each fusion junction candidate cluster: a. constructing a de-Bruijn graph with the k-mers of all cluster reads, b. executing bubble correction and tip removal algorithms, to correct sequencing errors, c. determining a number of contigs by traversing the de-Bruijn graph starting with a node that has the k-mer of the gene-specific primer that captured the fusion, and d. determining a longest contig as being a true assembled fusion transcript.

[0251] Here, the second predetermined set of metrics includes an assembled fusion transcript length, and the third predetermined set of metrics includes the number of unique transcript ends determined from the mapping output.

[0252] In this embodiment, categorizing further includes calculating a fusion score for each assembled fusion transcript based on any number of metrics of the first, second and third predetermined set of metrics, and the scoring function is calculated using the expression metrics, where m_ij is the value of the metric / for the assembled fusion transcript j, w_i is the weight associated with metric / , and n_\ is the expected mean of metric / , calculated from a control dataset with known fusions, comparing the fusion score to a first predetermined threshold and a second predetermined threshold, and categorizing the assembled fusion transcript as medium confidence if the fusion score is greater than the first predetermined threshold, categorizing the assembled fusion transcript as low confidence if the fusion score is less than the first predetermined threshold but greater than the second predetermined threshold, or categorizing the assembled fusion transcript as filtered if the fusion score is less than the second predetermined threshold, and outputting a result where low and medium confidence are subjected to orthogonal testing.

[0253] The metrics used for the scoring function include the read depth, the number of unique reads, the assembled fusion transcript length, and a ratio of read depth between fusion junction candidates and normal junction candidates with the same gene as the fusion target gene.

[0254] In a fourth embodiment of the first aspect, the predetermined number of base pairs to cluster junction candidates is 50 or less.

[0255] In a fifth embodiment of the first aspect, the method further includes determining a unique amplicon end by retrieving the last exon block mapping coordinate map_end and adding to it the length of the last mismatch block for each ofthe R2 reads. In this embodiment, identifying fusion junction candidates further includes identifying normal junction candidates in the merged mapping vector by searching for pairs of consecutive exon blocks from a same gene with a mismatch block in between, and recording the corresponding read pairs for each of the identified normal junction candidates. In this embodiment, the first predetermined set of metrics includes a number of unique amplicon ends, an estimated length of the fusion transcript, and a ratio between fusion junction candidates and normal junction candidates. The estimated length of the fusion transcript is determined by the average mapped length of all the read pairs in the corresponding junction candidate.

[0256] In a sixth embodiment of the first aspect, categorizing further includes categorizing any assembled fusion transcript as being an artifact if either the orientation of the fusion partners is not concordant or both partners are genes not targeted by the assay, or categorizing any assembled fusion transcript as being a supporting fusion if another assembled fusion with more reads is found on the same partner genes.

[0257] Example 2

[0258] An RNA-based amplicon assay for the sensitive and specific partner agnostic detection of fusions in poor quality samples involving 18 clinically relevant solid tumor driver genes is described. This assay performs well on fragmented FFPE RNA from clinical specimens with no prior knowledge of the fusion partner gene at low sequencing depths. Target enrichment is performed first to minimize loss of template information during adapter ligation. Then a random population of RNA splints act as scaffolds for annealing a universal adapter. After a second round of amplification and indexing, up to 21 clinical sample libraries are sequenced in multiplex. The actionable gene fusion content, cloud-based analysis, and integrative evidence algorithm provides accessible and focused results to ordering clinicians. These results are presented in a clinical reporting platform that enables rapid interpretation and reporting of clinical findings.

[0259] To address challenges associated with poor quality RNA, presequencing QC steps were employed. Imposing a qPCR CT threshold value of 30 for a housekeeping gene was found to be effective at screening samples for assay performance. Excluding low confidence fusion calls also improved specificity. With these filters, a sample rejection rate for archival clinical FFPE specimens with inputs as low as 10 ng of 32% (32 of 101) was observed. Despite having qPCR CT valuesabove 30, nine filtered samples still yielded correct fusion calls. Sample rejection was minimized by testing all samples and performing reflex testing with an alternative orthogonal DNA or antibody-based method on positive fusion calls in sub-threshold samples. This is advantageous due to the multi-target nature of NGS testing compared with immunohistochemical methods that are often single-target.

[0260] Successful performance was observed in archival FFPE specimens less than two years old down to 20 ng of input RNA, an amount obtainable from a single FFPE section. In an older cohort of up to five years of age, successful detection was observed with 100 ng of template. In this example, increasing template input beyond 100 ng did not improve detection of gene fusions in low quality samples.

[0261] Materials and Methods

[0262] Gene and exon target selection:

[0263] 18 driver genes (ALK, BRAF, CSF1, EGFR, ERG, FGFR1, FGFR2,FGFR3, MET, NRG1, NTRK1, NTRK2, NTRK3, PDGFB, PPARG, RAF1, RET, ROST) were targeted based on their presence and clinical actionability in solid tumors (Table 10).Table 10: Target genes associated with cancer type along with a non-limiting list of known / reported partner genes.

[0264] In total, 146 exons were selected based on biologically relevant configurations found in the literature and fusion databases1'3'26(Table 11).Table 11 : Target gene exon list. Target gene exons are listed for each reference transcript. Some exons described in the literature are present only in alternative transcripts and therefore are listed separately. The housekeeping gene TUBB is included for control analyses.

[0265] The tubulin gene TUBB was also included as a non-fusion housekeeping gene target to assess assay performance in all sample types.

[0266] Specimens and samples:

[0267] Universal normal human RNA (catalog R4234565-1 , Amsbio LLC, Abingdon, United Kingdom) was used as a fusion-negative analytical sample. The 18 fusion RNA mix (Seraseq catalog 0710-0497, LGC Clinical Diagnostics / SeraCare, Milford, Massachusetts, USA) was selected as a fusion-positive analytical sample. Double-stranded 500 bp geneblocks (gBIocks) (Integrated DNA Technologies (IDT), Coralville, Iowa, USA) were designed and centered around exonic junctions were also used as positive testing controls1 26. A 10 pg blend was spiked into 1 mg of the fusion-negative sample and used in runs (referred to as “74 Fusion-Positive”).Testing and validation specimens included FFPE cell pellets (catalog numbers: 3020- 0130 H2228, 3020-1130 H596, 3020-1430 HCC78, 3130-0130 U118MG, 3060-1130 A431 , Amsbio LLC, Abingdon, United Kingdom) and FFPE clinical samples sourced from commercial vendors with sample collection dates ranging from 2013-2019 (Reprocell, Beltsville, Maryland, Accio Biobank Online / Tissue for Research, United Kingdom, iSpecimen, Lexington, Massachusetts, USA, Trans-Hit Biomarkers Inc., Azenta, USA, Amsbio LLC, Abingdon, United Kingdom).

[0268] RNA Extraction and electrophoresis:

[0269] RNA was extracted from FFPE samples using a Maxwell RSC instrument and RNA FFPE kit (Promega, Madison, Wisconsin, USA). RNA yields were quantitated using the Qubit RNA kit (Invitrogen, Waltham, Massachusetts, USA). RNA DV200 was calculated using the Bioanalyzer RNA Pico assay (Agilent, Santa Clara, California, USA) as per manufacturer protocols.

[0270] cDNA Synthesis and quantitative PCR:

[0271] First strand cDNA was generated with random hexamers using the Superscript IV system (Invitrogen, Waltham, Massachusetts, USA). Second strand cDNA synthesis was generated using the New England BioLabs (NEB) ultra second strand synthesis module (NEB, Ipswich, Massachusetts, USA). The double stranded cDNA allows reverse or forward enrichment primers to be designed to enrich for 5’- partner-target-3’ or 5’-target-partner-3’ configurations. Using 1 ,8X AMPure XP beads (Beckman Coulter, Pasadena, California, USA), the purified cDNA was assessed for amplifiability with a PreSeq RNA QC assay as per manufacturer’s protocols (Archer Dx, In., Boulder, Colorado, USA).

[0272] Oligonucleotide design and dilution:

[0273] Nested enrichment primers and amplification primers were designed with Primer327. One enrichment and one amplification primer per exon was selected using a custom method based for multiplex optimization based on minimization of in silico free energy. Blocking oligos are the reverse-complement of enrichment primers. Enrichment primers, amplification primers, and blocking oligos were synthesized (Integrated DNA Technologies (IDT), Coralville, Iowa, USA) and pooled at 2 uM individual final concentration. The adapter and amplification primers contain 5’ Nextera XT extension sequences for creation of a Nextera XT sequencing library (Illumina, San Diego, California, USA).

[0274] Target Enrichment:

[0275] Double stranded cDNA was amplified with linear PCR using the Qiagen Multiplex PCR kit (Qiagen, Hilden, Germany) and a 2 uM multiplex enrichment primer mix in a 25 ul reaction. Thermal cycler program: 95°C 15 minutes, 8 cycles of (94°C 30 seconds, 60°C 3 minutes, 72°C 30 seconds), 4°C hold.Enrichment product is purified with 1 ,8X AMPure XP beads (Beckman Coulter, Pasadena, California, USA).

[0276] Adapter ligation:

[0277] Adapter and RNA splint pool were hybridized in annealing buffer (10mM Tris pH 8, 1 mM EDTA, 1 M NaCI) with thermal cycler program: 85°C 2 minutes, 80 cycles of 30 seconds with a reduction of -1 °C per cycle, 4°C hold. Enrichment products and blocking oligos were hybridized in annealing buffer with thermal cycler program: 95°C for 2 minutes, 70 cycles of 30 seconds with a reduction of -1 °C per cycle, 25°C hold. The hybridized duplexes were ligated using SplintR Ligase (New England BioLabs, Ipswich, Massachusetts, USA) with thermal cycler program: 25°C 8 hours, 9 cycles of (6 minutes, -1 °C per cycle), 16°C 30 minutes, 65°C 2 minutes, 4°C hold. Ligation product was purified with 1 ,8X AMPure XP beads (Beckman Coulter, Pasadena, California, USA).

[0278] Target Amplification:

[0279] The ligation product was amplified in KAPA HiFi hotstart readymix reaction (Roche, Basel, Switzerland) with a multiplex amplification primer pool, Nextera i7 R primer (targeting adapter) with thermal cycler program: 95°C 15 minutes, 8 cycles of (94°C 30 seconds, 60°C 3 minutes, 72°C 30 seconds), 4°C hold.Amplification product was purified with 1 ,8X Ampure XP beads (Beckman Coulter, Pasadena, California, USA).

[0280] Library Construction, Normalization and Sequencing:

[0281] A sequencing library was generated by PCR using the Nextera XT system (Illumina, San Diego, California, USA) using IDT for Illumina UDI adapters (Illumina, San Diego, California, USA) with KAPA HiFi PCR mix (Roche, Basel, Switzerland) using thermal cycler program: 95°C 30 seconds, 15 cycles (95°C 30 seconds, 55°C 30 seconds, 72°C 30 seconds), hold 4°C. Library product was purified with 0.7X AMPure XP beads(Beckman Coulter, Pasadena, California, USA). The libraries were quantified with a Qubit BR DNA kit (ThermoFisher Scientific, Waltham, Massachusetts, USA) and normalized to 20 nM. A pooled 4nM library is denatured, diluted, and sequenced with 5% PhiX using 2x150 cycle paired-end sequencing on a MiSeq™ instrument (Illumina, San Diego, California, USA).

[0282] Computational Methods:

[0283] Normal and cancer samples contain many chimeric reads that must be evaluated as biologically relevant fusions, artifacts, or alternative transcripts1 3. The sliding window k-mer based approach with two stage alignment for detection of atypical gene exons fusions paired with a classifier is executed on the input sequence data, which can be provided in FASTQ files. The following steps of the method correspond to the steps discussed earlier.

[0284] Read Filtering:

[0285] Reads are loaded from the FASTQ files, filtering out-of-target / mis- primed reads by labeling read pairs with primer identity based on the presence of the expected primer sequence at the start of the R1 read with maximum hamming distance of 1 from all the primer sequences in the assay manifest file. Labeled read pairs are aligned with the expected target sequence up to exon boundary. Alignments result in read pairs being assessed as on-target or off-target (when the primer binds to a different region) with tolerance for sequencing errors. Off-target reads are filtered and the remaining are used as an input for k-mer indexing.

[0286] K-mer indexing and read pseudo-alignment:

[0287] On-target read sequences are mapped to the reference human genome hg38 using a k-mer mapping algorithm. A k-mer index is created from on-target reads where a k-mer is a substring of length k of a given DNA sequence and is a structure where information about a given k-mer from a set of k-mers can be accessed (see Figure 5). The aim is to identify the corresponding exon regions of the sequencing read, to then find reads that map to two different genes (the driver gene in the known sequence, and a partner gene in the unknown sequence, see Figure 8), indicating a possible fusion supporting read. For a “perfect” fusion, where the entire sequence of the driver gene exon is immediately followed by the entire sequence of a exon from a different partner gene, we expect exactly k-1 k-mers covering the junction and with no discrepancies (gaps) between the mapped and expected exon lengths.

[0288] K-mer mapping and finding exon candidates

[0289] Given a k-mer from a sequencing read, a k-mer hit is the position of this k-mer in the read and the mapping information of this k-mer from the k-mer index. If a k-mer is not present in the index, the corresponding mapping information is called a miss. Repeating this process for all k-mers in a read results in a vector of k- mer hits, representing the original read. Consecutive hits for the same unique exon match are then aggregated. Consecutive misses are also merged. The combined misses are called a missed block. Exon candidates with a small range mapped are filtered, with a minimum mapped requirement of 24 bps. At the end of the process, a read is then represented as a mapping vector composed of consecutive exon candidates and misses.

[0290] After this step, additional steps of sequencing error / SNP correction and resolving mapping ambiguities as previously described are executed.

[0291] Unique amplicon determination:

[0292] This step is similar to the previously described step of amplicon end estimation. To estimate the number of unique amplicons on the list of supporting reads for all fusions, the aligned endpoints of the sequencing read R2 are used. Since the assay is based on an anchored linear amplification step, different sequenced amplicons will have different lengths, and therefore different endpoints (or starting R2 sequences). From the R2 read mapping, the end coordinate estimate is obtained by adding the last exon candidate mapping coordinate to the length of a missed region (if present, to account for errors at the end of the read).

[0293] Mate-pair Merging

[0294] After the k-mer mapping steps described above have been performed on both reads of a read pair, a merging algorithm combines exon candidates from both reads into a single mapping vector sorted by mapping coordinates. The end result is a single mapping vector, per read pair, used to find fusion events.

[0295] Junction Candidates

[0296] The merged mapping vectors are searched for junction candidates, defined as a pair of consecutive exon candidates on a given read mapping, with a miss region in between. If each exon is from a different gene (a known driver gene and an unknown partner gene), it will be called a fusion junction candidate, otherwise it is a normal junction candidate. Due to possible sequencing errors or SNVs close to the junction point, the previously described gap correction algorithm can be executed.

[0297] Clustering fusion junctions and recruiting supporting reads

[0298] For all junction candidates found in a read pair, if there is a single fusion junction candidate where the k-mer gap is smaller than a given threshold (default is k + 1), then the read pair will be added to a list of supporting reads for this fusion junction. If a read has no fusion junction candidates but it has normal junction candidates, it will be counted as a normal junction forthat particular gene. This is useful to estimate the ratio between the number of fusion supporting reads and normal reads for all fusion calls. After all read pairs are processed, all fusion junction candidates are clustered based on the genomic coordinates of the junction point, with a tolerance of 50 bps, meaning that junctions that are less than 50 bps apart are considered to be the same fusion.

[0299] Fusion Candidate Calling:

[0300] After the junction clustering step, for each fusion candidate, several metrics are calculated based on the supporting read pairs, such as the read depth, the number of unique amplicon ends, the estimated length of the fusion transcript, the ratio between fusion and normal reads, among other metrics.

[0301] Local assembly of fusion transcripts:

[0302] For each fusion, all its supporting reads are assembled using a de- Bruijn graph-based assembly. After the assembly, all reads of the sample aremapped against all assembled fusion transcripts to identify additional supporting reads that were filtered during the k-mer matching process.

[0303] Fusion scoring function

[0304] For a given fusion j, and a set of metrics, the fusion score s(j) is given by metrics

[0305] where m_ij is the value of the metric / for the fusion j, w_i is the weight associated with metric / , and p_\ is the expected mean of metric / , calculated from a control dataset with known fusions. The set of metrics selected were read depth, number of unique reads, fusion transcript length, and the fusion / normal read depth ratio.

[0306] Confidence Score Assignment

[0307] The fusion calling algorithm looks for reads with mapping regions in two different genes and must differentiate true fusions from a high number of artifact calls. All called fusions are scored based on a set of calculated metrics, and depending on the resulting score are either filtered, or receive a confidence label, from low, medium or high (see Figure 10). The fusion score is calculated using a set of metrics, weights, and thresholds calculated from a control dataset with known fusions with values defined as t_1 = 0 and t_2 = -1 (see Figure 11). The set of metrics selected were read depth, number of unique reads, fusion transcript length, and the fusion / normal read depth ratio. Fusions with score > t_1 are labeled medium confidence, while fusions with score < t_2 are filtered. The remaining fusions receive a low confidence label.1 . Artifact filter: The first step checks if the orientation of the fusion partners is concordant, and if one of the partners is a gene targeted by the assay. If either of the conditions is not satisfied, the fusion is filtered as an artifact.2. Supporting fusions: secondary fusion calls with the same exon junction between the same two genes are gathered as supporting of a primary fusion call with the highest support.3. Metrics filter: Several metrics are checked each with a minimum threshold value that the fusion has to pass. Failing on any of the metrics will label the call as filtered.4. High confidence: If a fusion call passed the previous steps and has read depth and number of unique reads higher than a given threshold, it is labeled high confidence.5. The remaining fusions are processed by a scoring function for reporting purposes. Depending on the score, a fusion is considered medium or low confidence, or be filtered.

[0308] Results

[0309] A method for enrichment and amplification of target gene exons

[0310] This example demonstrates a partner-agnostic NGS method and algorithm for targeted exon enrichment and detection of gene fusions in clinical FFPE tumors as shown in Figure 17, where the right side flowchart of Figure 17 “D: Fusion Analysis” is the same as Figure 4. Starting with FFPE tumor RNA template, the assay interrogated 146 exons in 18 driver genes for gene fusions. The assay employed random hexamer priming of cDNA to support fragmented and degraded FFPE RNA template. Double stranded cDNA was enriched for target exons using multiplex linear amplification from exon junction adjacent primers (Figure 18). The assay enabled target amplification without knowledge of the gene fusion partner by specifically ligating a universal adapter to the 3’ ends of enriched amplicons using a population of RNA splints with an adapter specific region and a random 10mer region along with the efficient Chlorella virus PBCV-1 DNA ligase28(Figure 19). Target amplification was performed using multiplex PCR using an adapter specific primer and an internally nested exon-specific amplification primer. Libraries were then generated, normalized, pooled, and sequenced on an Illumina MiSeq™ using paired- end 2 x 150 bps, yielding target exon sequence and junction coverage from the P5 read direction, and partner gene coverage from the P7 read direction (Figure 20).

[0311] One feature of the assay is the ability to detect unique template molecules based on the template cDNA molecule sequence since reads originating from the adapter have variable start site positions. To assess this, a normal human universal RNA template (fusion-negative) was used and the read coverage of the housekeeping gene TUBB (tubulin) was determined. TUBB exons 1-3 were targeted for amplification from within exon 4. Exons 1-3 had an average read depth of 4,279Xwhile the downstream regions (609-2,670) had an average read depth of 3X indicating effective and specific enrichment (Figure 21). From a sequencing run of 785,777 normal RNA derived reads, there were 12,588 TUBB gene transcript aligned reads, with 303 unique start sites.

[0312] Detection of gene fusions in reference samples

[0313] A custom analysis pipeline to identify on-target reads and the subset of candidate fusion reads was developed. On-target reads were collapsed into on- target unique reads based on adapter-read start site position relative to the aligned transcript. On-target reads with atypical exon organizations were classified as candidate fusion reads and scored. To establish specificity of read classification, a no-template control (NTC), a fusion-negative universal tissue normal human RNA control, and an 18 fusion positive control were assayed. The NTC contained no on- target reads, in contrast with the fusion-negative and the fusion-negative samples (Figure 22). On-target unique fusion reads were only identified in the 18 fusionpositive sample (median 12% from five runs), demonstrating enrichment of the initial ~3,000 copies per fusion target.

[0314] To establish the detection of fusions involving all 18 targeted genes 74 synthetic 500 bp dsDNA fusion templates (‘geneblocks’) were synthesized (Table 12).Table 12: Geneblock details. Called fusions all had junctions between end of upstream exon and start of downstream exon. Geneblocks not called had junctions between at least one exon starting within the exon or were reciprocal fusions that are not expected to be amplified.

[0315] Geneblocks encode biologically relevant fusions (65 two genes, 1 exon skipping) while 8 encoded non-biologically relevant fusions (2 reciprocal orientation, 6 mid-exon). 10 pg of the dsDNA geneblock pool spiked into the fusionnegative normal cDNA background were tested, and all 66 detectable geneblocks across 18 genes and 137 exons with the non-biologically relevant geneblock fusions filtered were identified. None of the geneblock fusions were detected in the fusionnegative samples.

[0316] Assessment of analytical sensitivity and specificity was performed using a dilution series experiment using the 18 fusion-positive sample. Three runs were performed representing ten replicates of the template mass input series: 12.5 ng, 25 ng, 50 ng, and 100 ng (750-12,000 copies per fusion). Unique reads supporting both normal and fusion transcripts increased with increasing input mass and fusion copies as provided by the manufacturer (Figure 23). Similar numbers of unique fusion reads were observed for each of the 18 fusion targets (Figure 24). Sensitivity of detection of fusions increased from 90.6% at 12.5ng (750 copies) to 98.3-99.4% for 25-100 ng (1 ,500-12,000 copies) with a minimum PPV of 99.4% (Table 13).Table 13: Sensitivity and PPV for the 18 fusion-positive sample with increasing input mass. True positives (TP), false negatives (FN), false positives (FP) are indicated. Sensitivity is defined as Sensitivity = (TP) / (TP+FN). PPV is defined as PPV = (TP) / (TP+FP).

[0317] Sensitivity of detection of SLC45A3-BRAF and EGFR-SEPT14 gene fusions in the contrived 18 fusion-positive sample was lowest (Table 14) with corresponding low fusion to normal read ratios (0.00-0.01 and 0.01-0.08) (Table 15).Table 14: Sensitivity = (TP) / (TP+FN) for expected fusions in the 18 fusionpositive sample with increasing template input mass. Assay of each input level was replicated for a total of 10 possible observations.Table 15: Median on-target read counts for target genes (includes any normal transcript or fusion transcript derived reads) for the NTC ‘no template control’, the fusion-negative sample, and 18 fusion-positive sample. Detected on-target fusion to normal read ratio range for observations is indicated for each fusion.

[0318] At 25 ng (~1500 copies per fusion) or higher, 16 of 18 fusions had100% detection. The fusions SLC45A3-BRAF and EGFR-SEPT14 were not detected in 6 of 60 observations at these input levels. These targets had high normal transcript reads and low fusion to normal read ratios of < 1 % and 1-8% respectively. Thefusions BRAF-SLC34A3, LMNA-NTRK1, TFG-NTRK1, and EGFR-SEPT14 were detected with 1-4% fusion to normal transcript read ratios but had lower numbers of normal transcripts. One false positive medium confidence artifact in the 18 fusionpositive sample dilution experiments for a single replicate (NFRKB-PDGFB) was observed. This fusion had two magnitudes lower read support compared to other fusions detected in this sample with 3 total and 1 unique supporting fusion read.

[0319] Detection of gene fusions in clinical FFPE specimens

[0320] To evaluate gene fusion detection from fusion-positive FFPE specimens of varying quality we performed a dilution series experiment on FFPE material from 5 treated cell lines (< 1 year old), and 15 archival clinical samples (7: < 2 years old, 8: 2-5 years old). RNA quality was assessed by DV200, and cDNA quality assessed using qPCR. RNA DV200 values ranged from 86-91% (median 90%) for the cell-lines and 30-86% (median 69%) for the archival specimens. Median qPCR CT values decreased with increasing template mass (20 ng, 100 ng, 250 ng, 500 ng), ranging from 26.6-23.3 CT value for cell lines and 33.2-27.8 CT value for archival specimens (Table 16).Table 16: Table of 20 FFPE testing samples and results.

[0321] Gene fusions were detected at all input levels for the FFPE cell lines and specimens less than 2 years old, while gene fusions were not detected at 20 ng inputs in the specimens that were older than 2 years (Figure 25). Increasing the input mass of RNA did not compensate for decreased quality.

[0322] Evaluation of sample quality with respect to assay performance was more effective with qPCR than DV200. At the lowest input of 20 ng of RNA, the sensitivity of detection was 65% without taking into account any quality assessment by DV200 and qPCR, 75% with DV200, and 100% with qPCR quality assessment (Table 17).Table 17: Sensitivity of the detection of 20 fusion-positive archival FFPE samples at different template input mass levels with and without screening by DV200 or qPCR.

[0323] A threshold CT value equal to 30 was used to reduce the inclusion of poor quality samples that can result in false negatives, however gene fusions were still detectable in some of these samples. There was no separation of data based on total unique reads and qPCR values with respect to template input amount showing that sample quality (qPCR) is more predictive of outcome than input mass (Figure 26).

[0324] The assay performance was next validated using 56 unique archival clinical FFPE samples that were represented across 101 observations (Table 18; four observations of a reference standard were omitted).Table 18: Table of validation samples and results. TN = True Negative, TP = True Positive, FP = False Positive, FN = False Negative.

[0325] Five validation sequencing experiments were performed by two different technologists on two instruments and on different days. Each sequencing run consisted of 3 controls and 21 samples. Samples were clinical FFPE specimens derived from patients with cancers of the lung, brain, lymph node, thyroid, and bone, with orthogonal molecular results including NGS testing, Immunohistochemistry (IHC), and Fluorescence In-situ Hybridization (FISH). When screened for minimum quality requirements (qPCR CT < 30), 69 observations from 44 unique samples were considered. There was no difference in the mass of template RNA used for low quality samples with qPCR CT > 30 (median 256 ng, minimum 11 ng) and high- quality samples with qPCR CT < 30 (median 223 ng, minimum 15 ng). Samples were classified as true positive, true negative, false positive, and false negative according to the assay relative to their orthogonal status.

[0326] Summary performance metrics were calculated with and without screening quality by qPCR as well as with and without low confidence fusion calls (Table19).Table 19: Assay validation performance metrics. Classifications at the sample level. Passing screen was defined as samples with qPCR CT values < 30 and / orwhere fusions were high or medium confidence as indicated. Accuracy of the assay where Accuracy = (TP+TN) / (TP+TN+FP+FN). Sensitivity of the assay where Sensitivity = (TP) / (TP+FN). Specificity of the assay where Specificity = (TN) / (TN+FP). Positive predictive value of the assay where PPV = (TP) / (TP+FP). Negative predictive value of the assay where NPV = (TN) / (TN+FN).

[0327] When all 101 sample observations from 56 samples were considered, accuracy was 85%, sensitivity was 80%, specificity was 97%, PPV was 98%, and NPV was 67%. Applying qPCR thresholds (CT < 30) decreased false negatives (improving NPV and sensitivity) while excluding low confidence fusion calls removed false positives (improving PPV and specificity). When both were considered, all metrics improved (accuracy = 94%, sensitivity = 90%, specificity = 100%, PPV = 100%, NPV = 87%).

[0328] There were four false negative observations after filtering by qPCR and fusion call confidence. One expected ROS1 fusion was not detected in a sample with low amplifiability near cut-off (qPCR CT 29.76) and had low unique read startsite diversity suggestive of low template quality. Samples with qPCR CT values below 30 typically have more than 300 unique read start-sites. Three additional samples had fusions called with low confidence and thus classified as fusionnegative. All three samples had low amplifiability (qPCR CT 28.16, 28.5, 28.5) and low fusion to normal read ratios for ALK (3%, 24%, 24%) which are typically high due to low normal transcript expression (average 92%, median 55% for all ALK fusion calls).

[0329] The predominant and most frequent unexpected fusion call observed in 5 clinical samples was SYN2-PPARG. These were observed in high quality (qPCR CT 23-25), high depth samples (5,221-7,746 on-target unique reads) and believed to represent a low level of read-through transcription between PPARG and adjacentgenes since false positive fusions involving the downstream genes TSEN2 and MKRN2 were also observed in other samples. Consequently, fusion calls between PPARG and SYN2, TSEN2, or MKRN2 are filtered. Detection of a SYN2-PPARG fusion has been reported in the literature without further investigation of clinical relevance29.

[0330] References

[0331] 1. Gao, Q. et al. Driver Fusions and Their Implications in theDevelopment and Treatment of Human Cancers. Cell Rep. 23, 227-238. e3 (2018).

[0332] 2. Yoshihara, K. et al. The landscape and therapeutic relevance of cancer-associated transcript fusions. Oncogene 34, 4845-4854 (2015).

[0333] 3. Jang, Y. E. et al. ChimerDB 4.0: an updated and expanded database of fusion genes. Nucleic Acids Res. 48, D817-D824 (2020).

[0334] 4. Stransky, N., Cerami, E., Schalm, S., Kim, J. L. & Lengauer, C.The landscape of kinase fusions in cancer. Nat. Commun. 5, 4846 (2014).

[0335] 5. Chambers, P. et al. Understanding Molecular Testing UptakeAcross Tumor Types in Eight Countries: Results From a Multinational Cross- Sectional Survey. JCO Oncol. Pract. 16, e770-e778 (2020).

[0336] 6. Mertens, F., Johansson, B., Fioretos, T. & Mitelman, F. The emerging complexity of gene fusions in cancer. Nat. Rev. Cancer 15, 371-381 (2015).

[0337] 7. Yun, J. W. et al. Dysregulation of cancer genes by recurrent intergenic fusions. Genome Biol. 21 , 166 (2020).

[0338] 8. Li, J. et al. A functional genomic approach to actionable gene fusions for precision oncology. Sci. Adv. 8, eabm2382 (2022).

[0339] 9. Du, Z. & Lovly, C. M. Mechanisms of receptor tyrosine kinase activation in cancer. Mol. Cancer 17, 58 (2018).

[0340] 10. Kumar-Sinha, C., Kalyana-Sundaram, S. & Chinnaiyan, A. M.Landscape of gene fusions in epithelial cancers: seq and ye shall find. Genome Med. 7, 129 (2015).

[0341] 11. Yamaoka, T., Kusumoto, S., Ando, K., Ohba, M. & Ohmori, T.Receptor Tyrosine Kinase-Targeted Cancer Therapy. Int. J. Mol. Sci. 19, E3491 (2018).

[0342] 12. Cohen, P., Cross, D. & Janne, P. A. Kinase drug discovery 20 years after imatinib: progress and future directions. Nat. Rev. Drug Discov. 20, 551 — 569 (2021).

[0343] 13. Schram, A. M., Chang, M. T., Jonsson, P. & Drilon, A. Fusions in solid tumours: diagnostic strategies, targeted therapy, and acquired resistance. Nat. Rev. Clin. Oncol. 14, 735-748 (2017).

[0344] 14. Laskin, J. et al. NRG1 fusion-driven tumors: biology, detection, and the therapeutic role of afatinib and other ErbB-targeting agents. Ann. Oncol. Off. J. Eur. Soc. Med. Oncol. 31 , 1693-1703 (2020).

[0345] 15. Ross, J. S. et al. The distribution of BRAF gene fusions in solid tumors and response to targeted therapy. Int. J. Cancer 138, 881-890 (2016).

[0346] 16. Zanwar, S. et al. Clinical and therapeutic implications of BRAF fusions in histiocytic disorders. Blood Cancer J. 12, 97 (2022).

[0347] 17. Nikanjam, M., Okamura, R., Barkauskas, D. A. & Kurzrock, R.Targeting fusions for improved outcomes in oncology treatment. Cancer 126, 1315— 1321 (2020).

[0348] 18. Weiss, L. M. & Funari, V. A. NTRK fusions and Trk proteins: what are they and how to test forthem. Hum. Pathol. 112, 59-69 (2021).

[0349] 19. Schmitt, F., Di Lorito, A. & Vielh, P. Molecular Testing onCytology for Gene Fusion Detection. Front. Med. 8, 643113 (2021).

[0350] 20. Tachon, G. et al. Targeted RNA-sequencing assays: a step forward compared to FISH and IHC techniques? Cancer Med. 8, 7556-7566 (2019).

[0351] 21 . Wu, Y.-C. et al. Comparison of IHC, FISH and RT-PCRMethods for Detection of ALK Rearrangements in 312 Non-Small Cell Lung Cancer Patients in Taiwan. PLOS ONE 8, e70839 (2013).

[0352] 22. Somaschini, A. et al. Mining potentially actionable kinase gene fusions in cancer cell lines with the KuNG FU database. Sci. Data 7, 420 (2020).

[0353] 23. Bruno, R. & Fontanini, G. Next Generation Sequencing forGene Fusion Analysis in Lung Cancer: A Literature Review. Diagnostics 10, 521 (2020).

[0354] 24. Heydt, C. et al. Detection of gene fusions using targeted nextgeneration sequencing: a comparative evaluation. BMC Med. Genomics 14, 62 (2021).

[0355] 25. Schroder, J., Kumar, A. & Wong, S. Q. Overview of FusionDetection Strategies Using Next-Generation Sequencing. Methods Mol. Biol. Clifton A / J 1908, 125-138 (2019).

[0356] 26. Tate, J. G. et al. COSMIC: the Catalogue Of SomaticMutations In Cancer. Nucleic Acids Res. 47, D941-D947 (2019).

[0357] 27. Untergasser, A. et al. Primer3-new capabilities and interfaces.Nucleic Acids Res. 40, e115 (2012).

[0358] 28. Lohman, G. J. S., Zhang, Y., Zhelkovsky, A. M., Cantor, E. J. &Evans, T. C. Efficient DNA ligation in DNA-RNA hybrid helices by Chlorella virus DNA ligase. Nucleic Acids Res. 42, 1831-1844 (2014).

[0359] 29. Soon, G. S. T., Chang, K. T. E., Kuick, C. H. & Petersson, F. A case of nasal low-grade non-intestinal-type adenocarcinoma with aberrant CDX2expression and a novel SYN2-PPARG gene fusion in a 13-year-old girl. Virchows Arch. Int. J. Pathol. 474, 619-623 (2019).

[0360] 30. Latysheva, N. S. & Babu, M. M. Discovering and understanding oncogenic gene fusions through data intensive computational approaches. Nucleic Acids Res. 44, 4487-4503 (2016).

[0361] 31 . Groelz, D., Viertler, C., Pabst, D., Dettmann, N. & Zatloukal, K.Impact of storage conditions on the quality of nucleic acids in paraffin embedded tissues. PLoS ONE 13, e0203608 (2018).

[0362] 32. Malnar, M. & Rezen, T. Factors affecting RNA quantification from tissue long-term stored in formalin. J. Pharmacol. Toxicol. Methods 96, 61-66 (2019).

[0363] 33. von Ahlfen, S., Missel, A., Bendrat, K. & Schlumpberger, M.Determinants of RNA Quality from FFPE Samples. PLoS ONE 2, e1261 (2007).

[0364] 34. Kresse, S. H. et al. Evaluation of commercial DNA and RNA extraction methods for high-throughput sequencing of FFPE samples. PLoS ONE 13, e0197456 (2018).

[0365] 35. Cieslik, M. et al. The use of exome capture RNA-seq for highly degraded RNA with application to clinical cancer sequencing. Genome Res. 25, 1372-1381 (2015).

[0366] All citations are hereby incorporated by reference.

[0367] The present invention has been described with regard to one or more embodiments. However, it will be apparent to persons skilled in the art that a number of variations and modifications can be made without departing from the scope of the invention as defined in the claims. Therefore, although various embodiments of the invention are disclosed herein, many adaptations and modifications may be made within the scope of the invention in accordance with the common general knowledge of those skilled in this art. Such modifications include the substitution of known equivalents for any aspect of the invention in order to achieve the same result in substantially the same way. Numeric ranges are inclusive of the numbers defining the range. In the specification, the word “comprising” is used as an open-ended term, substantially equivalent to the phrase “including, but not limited to,” and the word “comprises” has a corresponding meaning. It is to be however understood that, where the words “comprising” or “comprises,” or a variation having the same root, are used herein, variation or modification to “consisting” or “consists,” which excludes any element, step, or ingredient not specified, or to “consisting essentially of’ or “consists essentially of,” which limits to the specified materials or recited steps together with those that do not materially affect the basic and novel characteristics of the claimed invention, is also contemplated. Citation of references herein shall not beconstrued as an admission that such references are prior art to the present invention. All publications are incorporated herein by reference as if each individual publication was specifically and individually indicated to be incorporated by reference herein and as though fully set forth herein. The invention includes all embodiments and variations substantially as hereinbefore described and with reference to the examples and drawings.

Claims

WHAT IS CLAIMED IS:1 . A method for ligating single stranded DNA molecules comprising: i) providing a first single stranded DNA molecule comprising a known sequence of at least 10 nucleotides in length; ii) providing a second single stranded DNA molecule comprising an unknown sequence of least 10 nucleotides in length; iii) providing a single stranded RNA molecule comprising a 5’ end substantially complementary to the known sequence of the first single stranded DNA molecule and a 3’ end comprising a random sequence of least 10 ribonucleotides in length; iv) providing a DNA ligase; and v) combining the first single stranded DNA molecule, the second single stranded DNA molecule, the single stranded RNA molecule and the DNA ligase, wherein the first single stranded DNA molecule and the second single stranded DNA molecule are ligated to provide a ligated single stranded nucleic acid molecule when the unknown sequence of the first single stranded DNA molecule is substantially complementary to the random sequence of the single stranded RNA molecule.

2. The method of claim 1 wherein the first single stranded DNA molecule is blocked at the 3’ end.

3. The method of claim 1 or 2 wherein the first single stranded DNA molecule further comprises a spacer sequence adjacent to the 3’ end of the known sequence.

4. The method of any one of claims 1 to 3 wherein the known sequence of the first single stranded DNA molecule is from ten to 100 nucleotides in length.

5. The method of any one of claims 1 to 4 wherein the unknown sequence of the second single stranded DNA molecule is from ten to 10,000 nucleotides in length.

6. The method of any one of claims 1 to 5 wherein the single stranded RNA splint molecule is about 20 nucleotides in length and wherein 10 nucleotides at the 5’ end are substantially complementary to the known sequence of the first singlestranded DNA molecule and the random sequence is substantially complementary to the unknown sequence of second single stranded DNA molecule.

7. The method of any one of claims 1 to 6 wherein the first single stranded DNA molecule, the second single stranded DNA molecule, and the single stranded RNA molecule are combined under conditions suitable for annealing of the first single stranded DNA molecule and the second single stranded DNA molecule to the single stranded RNA molecule.

8. The method of any one of claims 1 to 7 further comprising an antiprimer molecule, wherein the first single stranded DNA molecule is combined with the single stranded RNA molecule and the second single stranded DNA molecule is combined with the antiprimer sequence prior to combination of the first single stranded DNA molecule and the single stranded RNA molecule duplex and the second single stranded DNA molecule and antiprimer duplex.

9. The method of any one of claims 1 to 8 wherein the DNA ligase is a Chlorella virus DNA ligase.

10. The method of any one of claims 1 to 9 wherein the second single stranded DNA molecule is in vitro ssDNA or single stranded cDNA.11 . The method of any one of claims 1 to 10 wherein the second single stranded DNA molecule further comprises a second known sequence 5’ to the unknown sequence.

12. The method of any one of claims 1 to 11 wherein the second single stranded DNA molecule further comprises a primer sequence at the 5’ end.

13. The method of claim 12, further comprising amplifying the ligated single stranded nucleic acid molecule.

14. The method of claim 13 wherein the amplifying is performed by PCR.

15. The method of claim 14 wherein the PCR is linear PCR or nested PCR.

16. A kit comprising: i) a first single stranded DNA molecule comprising a known sequence of at least 10 nucleotides in length;ii) a second single stranded DNA molecule comprising an unknown sequence of least 10 nucleotides in length; and iii) a single stranded RNA splint molecule comprising a 5’ end substantially complementary to the known sequence of the first single stranded DNA molecule and a 3’ end comprising a random sequence of least 10 ribonucleotides in length.

17. The kit of claim 16 further comprising a DNA ligase.

18. The kit of claim 17 wherein the DNA ligase is a splint ligase.

19. The kit of any one of claims 16 to 18 further comprising instructions for use.

20. A computer implemented method for identifying from a plurality of read pairs R1 and R2 of a sequenced sample with inserts derived from a set of PCR amplified adapter ligated DNA molecules, any target gene fused with any other DNA sequence which are potentially clinically actionable, comprising: generating a k-mer index of predetermined length k from a reference genome or a transcriptome; mapping each sequencing read from the plurality of read pairs to exons of the reference genome with the k-mer index to identify corresponding exon regions of the sequencing read, and generating for each sequencing read a mapping vector including one or more blocks, where each block is either an exon block defined by a range of consecutive k-mer matches and their corresponding locations on an exon of the reference genome, or a mismatch block defined by a length of consecutive k-mer mis-matches; identifying exon blocks with overlapping consecutive k-mer matches in each of the mapping vectors, and merging the same exon blocks with adjustment of mapping coordinates or retaining the exon block having predetermined high priority criteria to resolve mapping ambiguities; combining both the mapping vectors corresponding to each read pair and sorting by genomic coordinates to provide intermediate mapping vectors, and merging any exon blocks with the same respective exons with adjustment of mapping115coordinates to generate merged mapping vectors; identifying fusion junction candidates in the merged mapping vectors by searching for a pair of consecutive exon blocks from different genes, called fusion partner genes, with a mismatch block in between, and recording the corresponding read pairs for each of the identified fusion junction candidates; clustering all the fusion junction candidates based on genomic coordinates of junction points of the pairs of consecutive exons from different genes where junctions closer than a predetermined number of base pairs are considered to be the same fusion junction candidate, and calculating a first predetermined set of metrics based on the read pairs including read depth; assembling a fusion transcript for each fusion junction candidate cluster using the corresponding recorded read pairs to generate assembled transcripts, and calculating a second predetermined set of metrics from the assembled transcripts; mapping the plurality of read pairs to the fusion transcripts, and calculating a third predetermined set of metrics based on the mapping output including a number of mapped reads; and categorizing any fusion transcript as being a high confidence reportable fusion when the orientation of both fusion partners is concordant, exactly one of the partners is a gene targeted by the assay, minimum thresholds for specific metrics of the first, second and third predetermined set of metrics are met, and the number of mapped reads and the read depth exceed predetermined thresholds.21 . The method of claim 20, wherein exon blocks having less than a minimum mapped bps are filtered out.

22. The method of claim 21 , wherein the minimum mapped bps is 24 bps.

23. The method of claim 20, wherein each of the mapping vectors is represented116by the expression (e_n, read_st— >read_end, map_st^map_end) for each exon block and (miss, miss_st— >miss_end, len) for each mismatch block, where e_n is the exon unique identifier, read_st and miss_st are the starting base pair on the sequencing read, read_end and miss_end are the ending base pair on the sequencing read, map_st is the starting mapped region on the exon for read_st, map_end is the ending mapped region on the exon for read_end, miss is a status identifier for a mis-match, and len is the bps length of the mismatch block.

24. The method of claim 20, further including error correcting the mapping vectors to remove mismatch blocks interposed between two exon blocks with the same exon.

25. The method of claim 24 wherein error correcting includes removing the length of consecutive k-mer mis-matches of the mismatch block interposed between a first range of consecutive k-mer matches of a first exon block and a second range of consecutive k-mer matches of a second exon block of the same exon, when a continuation of the first range by the length of consecutive k-mer is sequentially followed by a beginning of the second range.

26. The method of claim 25, wherein error correcting to remove a mismatch block between two consecutive exon blocks with the same exon includes determining when read_end + len + 1 for the first exon block = read_st for the subsequent exon block to remove the mismatch block.

27. The method of claim 25, wherein error correcting to remove the length of consecutive k-mer mis-matches includes determining when map_end + len + 1 for a first exon = map_st for a subsequent exon to remove the length of consecutive k-mer mis-matches.

28. The method of claim 23, wherein identifying fusion junction candidates includes executing an error correction applied on k-mer gaps g_1 , g_2 and g_k of the fusion junction candidates when g_k > 0, where g_1 is a distance between map_end of the first exon block to the corresponding exon extremity, g_2 is a distance between map_st of the second exon block to the corresponding exon extremity, g_k is the difference between the k-mer length k and the mismatch block len, and the k-mer gap correction steps include while g_k > 0117a. if g_1 > 0 then i. g_k = g_k - 1H. g_i = g - 1 b. else If g_2 > 0 then i. g_k = g_k -1 ii. g_2 = g_2 - 129. The method of claim 23, wherein the predetermined high priority criteria include an exon block with a larger mapping range than the other, an exon block where its exon comes from an assay targeted gene, an exon block matching a list of preferred candidates, or when both exons are from the same gene, selecting the exon with highest priority, calculated from sorting the corresponding reference genome transcripts for each exon in descending order of length.

30. The method of claim 29, wherein if none of the predetermined high priority criteria is met, then remove the exon block having the largest exon gap, where an exon gap is a distance between map_end of a first exon and map_st of a subsequent exon to their respective exon extremities.31 . The method of claim 20, wherein the predetermined number of base pairs to cluster junction candidates is 50 or less.

32. The method of claim 20, further including determining a unique amplicon end by retrieving the last exon block mapping coordinate map_end and adding to it the length of the last mismatch block for each of the R2 reads.

33. The method of claim 32, wherein identifying fusion junction candidates further includes identifying normal junction candidates in the merged mapping vector by searching for pairs of consecutive exon blocks from a same gene with a mismatch block in between, and recording the corresponding read pairs for each of the identified normal junction candidates.

34. The method of claim 33, wherein the first predetermined set of metrics includes a number of unique amplicon ends, an estimated length of the fusion118transcript, and a ratio between fusion junction candidates and normal junction candidates.

35. The method of claim 34, wherein the estimated length of the fusion transcript is determined by the average mapped length of all the read pairs in the corresponding junction candidate.

36. The method of claim 30, wherein assembling the fusion transcript includes the following steps for each fusion junction candidate cluster a. constructing a de-Bruijn graph with the k-mers of all cluster reads, b. executing bubble correction and tip removal algorithms, to correct sequencing errors, c. determining a number of contigs by traversing the de-Bruijn graph starting with a node that has the k-mer of the gene-specific primer that captured the fusion, and d. determining a longest contig as being a true assembled fusion transcript.

37. The method of claim 36, wherein the second predetermined set of metrics includes an assembled fusion transcript length.

38. The method of claim 37, wherein, the third predetermined set of metrics includes the number of unique transcript ends determined from the mapping output.

39. The method of claim 20, wherein categorizing further includes categorizing any assembled fusion transcript as being an artifact if either the orientation of the fusion partners is not concordant or both partners are genes not targeted by the assay.

40. The method of claim 20, wherein categorizing further includes categorizing any assembled fusion transcript as being filtered if minimum thresholds for any one of the specific metrics is not met.41 . The method of claim 20, wherein categorizing further includes categorizing any assembled fusion transcript as being a supporting fusion if another assembled fusion with more reads is found on the same partner genes.

42. The method of claim 38, wherein categorizing further includes119calculating a fusion score for each assembled fusion transcript based on any number of metrics of the first, second and third predetermined set of metrics, and the scoring function is calculated using the expression , tor i t metricsi , where m_ij is the value of the metric / for the assembled fusion transcript j, w_i is the weight associated with metric / , and n_\ is the expected mean of metric / , calculated from a control dataset with known fusions, comparing the fusion score to a first predetermined threshold and a second predetermined threshold, and categorizing the assembled fusion transcript as medium confidence if the fusion score is greater than the first predetermined threshold, categorizing the assembled fusion transcript as low confidence if the fusion score is less than the first predetermined threshold but greater than the second predetermined threshold, or categorizing the assembled fusion transcript as filtered if the fusion score is less than the second predetermined threshold, and outputting a result where low and medium confidence are subjected to orthogonal testing.

43. The method of claim 42, wherein the metrics used for the scoring function include the read depth, the number of unique reads, the assembled fusion transcript length, and a ratio of read depth between fusion junction candidates and normal junction candidates with the same gene as the fusion target gene.120

Images