A single CML case from a patient of the Jiang Yin People’s Hospital forms the basis of this report. The sample for gene fusion study was obtained in March 2011 and the patient was a Chinese male aged 78 years old. Clinical information was obtained from patient charts in the Jiang Yin People’s Hospital. CML in CP was diagnosed in 2005 on the basis of a morphological study of bone marrow and peripheral blood specimens (total of 8.1% myeloblast and progranulocyte in blood smear with WBC=34.8×10⁹/l, Hb=62 g/l and PLT=529×10⁹/l). In addition to increased myeloblast and progranulocyte levels, additional symptoms were identified, including blood dilution and splenomegaly demonstrated using ultrasonography. However, b3a2 b2a2 and ela2 fusion junctions in the BCR-ABL gene were not detected in 2008. One year later, detection of the classic BCR-ABL breakpoint remained negative.

Written informed consent from the patient was obtained and the present study was reviewed and approved by the Clinical Laboratories Department of the People’s Hospital in Jiang Yin City (Jiang Yin, China).

Total RNA was extracted from the peripheral blood using TRIzol according to the manufacturer’s instructions (Invitrogen, Carlsbad, CA, USA). Following the TruSeq RNA SamplePrep Guide (Illumina, Inc., San Diego, CA, USA), the Illumina standard kit was used for mRNA-Seq sample preparation. Briefly, 10 μg of total RNA extracted from the sample was used for polyA mRNA selection with polyT oligo-conjugated magnetic beads by two rounds of purification, followed by thermal mRNA fragmentation. Using reverse transcriptase (SuperScript II) and random primers, the fragmented mRNA was subjected to cDNA synthesis. Following this, cDNA was converted into double-stranded cDNA and after end repair [Klenow fragment, T4 polynucleotide kinase, T4 polymerase and 3-‘A’ add process (Klenow exo-fragment)]the product was ligated to Illumina Truseq adaptors. A 2% agarose gel was used to perform size selection, which generated 380-bp cDNA libraries. Finally, the libraries were enriched through 15 cycles of PCR and purified with the QIAquick PCR purification kit (Qiagen, Hilden, Germany). Elution buffer was used to dilute the enriched libraries to a final concentration of 10 nM.

Libraries from CML blood in CP were analyzed at a concentration of 11 pM on a single Genome Analyzer IIx (GAIIx) lane using 115-bp sequencing. Raw RNA-Seq data were filtered using FASTx-tools (http://hannonlab.cshl.edu/fastx_toolkit/) according to the following criteria: i) reads containing sequencing adaptors were removed; ii) nucleotides with a quality score <20 were trimmed from the end of the sequence; iii) reads <50 were discarded; and iv) artificial reads were removed. Following the filtering pipeline, all clean reads were trimmed to 90 bp long due to 3′ low quality base. A total of 3.2 Gbp (3,193,208,820) of cleaned, paired-end reads was produced. Raw sequence data were submitted to the NCBI Short Read Archive (accession number, SRP011486).

Clean and trimmed reads were aligned with the Ensembl H. sapiens reference genome (build GRCh37) using TopHat v1.3.3 (14).The program initially removes a portion of the reads based on quality information accompanying each read and then aligns the reads to the reference genome. TopHat allows multiple alignments per read (up to 20 by default) and a maximum of two mismatches when mapping the reads to the reference genome. The default parameters were used. Aligned read files were then processed by Cufflinks v1.2.1 (15), which measures the relative expression of the genes with the normalized RNA-Seq fragment counts. Fragments per kilobase of exon per million fragments mapped (FPKM) is used as the unit of measurement. A Bayesian inference method (16) was used to calculate confidence intervals for FPKM estimates. The reference GTF annotation file used in Cufflinks was downloaded from the Ensembl database (Homo_sapiens.GRCh37.63.gtf) (17). The gene expression data were submitted to the GEO database (accession ID, GSE36522).

Using TopHat, all filtered RNA-Seq reads were mapped to the reference transcriptome, which was downloaded from the Ensembl database (Homo_sapiens.GRCh37.63.cdna.all.fa). Read pairs mapping to the same transcripts were removed. Remaining reads were analyzed by two software packages, deFuse (deFuse-0.4.2) (18) and TopHat-Fusion (TopHat-Fusion-0.1.0) (19), to identify candidate gene fusions. The bowtie-index used in the TopHat-Fusion was downloaded from the TopHat homepage (H. sapiens UCSC hg19). The parameters of the TopHat-Fusion used were obtained from the ‘Getting Started’ (http://tophat-fusion.sourceforge.net/tutorial.html) tutorial. The deFuse parameters used were default.

deFuse results were filtered according to McPherson et al(18) and Steidl et al(20). Two candidate gene fusions were obtained following the filtering pipeline. The TopHat-Fusion results were parsed by in-house Perl scripts and filtered according to the following pipeline: i) span reads were >8 reads; ii) ratio of against reads vs. span reads was <0.5; and iii) gene fusions involving ribosomal proteins or small nuclear ribosomal proteins were excluded. Following parsing, only one filtered candidate fusion remained. Filtered candidates detected simultaneously by deFuse and TopHat-Fusion were considered reliable candidate gene fusions. Following this filtering pipeline, a single reliable candidate, RNF213-SLC26A11, was obtained.

To reduce the false-positive rate and remove read-throughs, we improved the filtering standard by comparing the conditions of our filtering method with those of the final results. Two filtering standards of deFuse were modified: i) the number of splitr_count was revised according to the average coverage of the sequencing depth and ii) the events, which happened between 2 adjacent genes without any eversion or inversion, were considered to be produced by read-through. The detailed conditions of the 14 candidates are listed in Table I.

To detect fusion transcripts, we designed a forward primer targeting the 5′ partner gene and reverse primer targeting the 3′ partner. Primer pairs (Table II) for the coding exons of the fusion genes were generated using Primer 5 software (Premier Biosoft International, Palo Alto, CA, USA) and the PCR volume comprised 20 μl sample, 2 μl 10X PCR buffer, 2 μl cDNA template, 0.4 μl dNTP, 0.4 μl Taq enzyme (Genscript, Piscataway, NJ, USA) and 0.2 pmol/μl each oligonucleotide. PCR was performed using the following procedure: 95°C for 2 min, 35 cycles of 95°C for 15 sec, 62°C for 20 sec and 72°C for 20 sec, followed by 72°C for 2 min. The presence of the fusion gene in the CML blood sample was confirmed. PCR products of the fusion gene were cloned in the pGEM^®-T Easy Vector (Promega Corp., Madison, WI, USA) and then sequenced with the T7 primer using an Applied Biosystems 3730 DNA Analyzer (Life Technologies Corporation, USA).

Views of the candidate gene fusions were generated by running CIRCOS (available at http://circos.ca/). Reads coverage of the candidate gene fusions in the sample were observed in the Integrative Genomics Viewer (IGV).

CML blood in CP was obtained from a male patient of 78 years old. The peripheral blood smear is demonstrated in Fig. 1. Massively parallel paired end cDNA sequencing was performed using the Illumina Genome Analyzer IIx. In total, 43.6 million raw reads were produced. Following filtering and trimming, 35.5 million clean reads were obtained. All reads were aligned to the reference genome by TopHat. In summary, 15.3 million read pairs were uniquely aligned. Reads (74%) were mapped to Ensembl reference genes. The average coverage of the sequencing depth was >30-fold of the human transcriptome. Only 4% reads were mapped to the H. sapiens mitochondrial genome. Raw 115-bp long reads from a single lane on the platform of Genome Analyzer IIx (GAIIx) was analysed. In addition, 20.89% (7,413,014) reads were splice junction reads. All observations were indicative of a correctly constructed library.

To measure gene expression in the sample, we utilized Cufflinks to estimate gene expression. Normalized expression levels of each gene were measured in FPKM. By selecting genes with an FPKM value >1, 22,607 expressed genes were detected in the sample, including the majority of annotated reference genes in H. sapiens.

Two algorithms, including deFuse and TopHat-Fusion, were used to detect gene fusion events in the sample. Following removal of read pairs mapping to the same genes, 5.4 million (~30.9%, 5.4/17.7) reads pairs were obtained for gene fusion event prediction. Several gene fusion candidates were discovered, the majority of which were located within the chromosome. Following a manual check of the raw results produced by deFuse, unreliable results were removed. The filtered results of TopHat-Fusion and deFuse were combined to produce 14 gene fusion candidates (Table I). Among the candidates, one gene fusion event RNF213-SLC26A11 was detected by deFuse and TopHat-Fusion with a 3-bp breakpoint shift.

The majority of candidate gene fusions were identified within the same chromosome (Table I). Using manual checks and selection, we identified that the majority of fusion events occurred between adjacent genes on the same chromosome, which were usually regarded as read-through (21), while only two candidates, including RNF213-SLC26A11 and BCR-ABL1, were the most reliable gene fusion pairs. BCR-ABL is a well-documented gene fusion in CML and is used as a standard method for CML diagnosis (22). The RNF213-SLC26A11 fusion has not been previously identified in CML. In the present study, this novel fusion was further validated by reverse transcription polymerase chain reaction (RT-PCR) and the function of RNF213-SLC26A11 was predicted using bioinformatic methods.

To experimentally confirm the novel gene fusion identified by RNA-Seq, the expression level of RNF213-SLC26A11 in the sample was validated by RT-PCR. A primer pair, located at the second intron of RNF213 (5′-GACTCCTGCTCTTGCTTC TGG-3′) and the 8th exon region of SLC26A11 (5′-ATCGTC CCGTTGGCTGTG-3′) was designed. The results confirmed the existence of the fusion event in the sample (Fig. 2C), consistent with conclusions based on RNA-Seq analysis.

As demonstrated in Fig. 2A, RNF213 and SLC26A11 are adjacent gene pairs located in chromosome 17 separated by 7 kbp. The genes are transcribed in the same direction in wild-type individuals. In the present sample, we observed a fusion between the second intron of RNF213 and the eighth exon of SLC26A11, generating a chimeric RNF213-SLC26A11 transcript, which has not previously been identified.

Analysis of reads coverage of RNF213 and SLC26A11 demonstrated that each exon of RNF213 was expressed at normal levels and the first 7 exons of SLC26A11 were expressed at extremely low levels (Fig. 2B), indicating that i) the RNF213-SLC26A11 fusion may occur in the majority of cells or ii) the fusion event altered expression of SLC26A11 using the promoter of RNF213, largely caused by the duplication of the first 2 exons of RNF213. RNF213 was normally expressed. In addition, compared with domains of the whole SLC26A11 gene under normal conditions (Fig. 3A), specific domains of the gene SLC26A11 in the chimeric transcript, including Sulfate_tra_GLY (pfam13792), EF_0839 (TIGR03581), sulfate_transp (pfam00916), sulP (TIGR00815), SUL1 (COG0659) and PRK11660 (PRK11660), were partially damaged (Fig. 3B). These observations indicate that the normal function of RNF213 did not appear to be affected, while loss of function in SLC26A11 was identified. However, further investigation is required to understand the specific mechanism and functional consequence of this fusion event.

In addition to the novel gene fusion, we identified another gene fusion candidate, which had been previously identified in CML (23). As demonstrated in Fig. 4B, BCR is located in chromosome 22 while ABL1 is on chromosome 9. In the present sample, a fusion was identified between the first intron of BCR and the second exon of ABL1, generating a chimeric transcript BCR-ABL1 (Fig. 4). Deininger et al(24) previously reported that this fusion is a e1a2 fusion junction and is extremely rare (25).