PVT: An Efficient Computational Procedure to Speed up Next-generation Sequence Analysis

TopHat and its pitfalls

The different steps along with their corresponding functions involved in the ‘spliced alignment’ is discussed in details in Methods section (summarised in Additional file 3: Table S3). The order of execution of these steps during ‘spliced alignment’ using TopHat is shown in Additional file 4: Figure S1. The execution time for the single end read input dataset using TopHat is shown in Additional file 5: Figure S2. Computational resource (i.e. CPU and memory) utilization at different steps of ‘spliced alignment’ was noted for the entire single end read dataset. We have shown the result for SRR094775 data only in Figure 1, since this data consists of maximum number of filtered reads.

For the find_juncs step, we found that the CPU remains significantly underutilized (~12% utilization) for a considerable span of time as compared to other steps. Moreover, the average CPU utilization throughout the entire execution period for a single dataset is ~39%. Further cache memory utilized by TopHat during find_juncs step is also more compared to other steps These observations from the utilization graphs (Figure 1) gave us the clue that find_juncs is that step which can be modified to bring significant improvement in the efficiency of TopHat. Hence, parallelizing find_juncs would lead to an increase in CPU and memory utilization and decrease execution time. However, care should be taken towards the extent of parallelization so as to ensure no buffer overflow or CPU overutilization.

While analyzing the steps and workflow of ‘spliced alignment’, in details we found a redundant step in TopHat where the alignment tool builds the bowtie indices for aligning the reads with the genome repeatedly for processing each data file of the same experiment. This adds on to the resource consumption. These indices can be built once and used for different data files unless the reference genome is updated or a new reference organism is being worked upon. We could increase the efficiency of TopHat and reduce its execution time by building the bowtie reference indices only once for the analysis of the entire data set of an experiment.

PVT

The performance analysis reports obtained on executing TopHat for single end read dataset (as shown in Figure 1) motivated us to develop PVT (Pipelined Version of TopHat) wherein we modified the execution workflow of TopHat for efficient utilization of CPU, memory and time. In PVT, the steps filter_reads, gene_align and genome_align are executed serially but we parallelized the find_juncs step so as to ensure proper utilization of computational resources, thereby reducing the execution time. Figure 2 shows the PVT alignment pipeline and its order of execution for single end reads. Rather than searching for the plausible junctions in the whole reference genome at a go, we looked for the junctions simultaneously across all reference chromosomes of the genome. Thus we were able to parallelize find_juncs by first splitting the reads that mapped to the genome (i.e. the output of genome_align), based on chromosomal reference and then concatenating the split files into comparable chromosomal reference (using splitBy_chromReference). The find_juncs is then run on separate chromosomal file references in parallel and the output of possible junctions are concatenated (using concatenate_segments). This is followed by serial execution of the steps- junc_align, span_reads and report. Parallel execution was achieved by executing the steps in both foreground and background. The blue and green bordered boxes in Figure 2 represent the processes run in the foreground and background respectively. The parallelization of this step and the execution of the subsequent steps is scheduled in such a way that no dependent step begins its execution before the completion of the previous step.

Prior to executing an independent step of the pipeline, clearing cache facilitates further parallelization. In PVT, the cache memory is cleared at every step to save the memory and meet the demand of other steps with high memory requirement. This helps avoid thrashing due to memory shortage and hence brings down the memory requirement for find_juncs (as shown in Figure 3).

In PVT, we have also addressed the problem of additional resource consumption for building bowtie indices. Here, the bowtie indices were built just once for analyzing the entire dataset, thereby removing the redundancy involved in the process followed by TopHat. This indexing took ~10 minutes in our standalone configuration for the single end read analysis, which reduces 6.25% of the total execution time taken by TopHat.

Improvement in find_juncs: The CPU utilization during the execution of PVT pipeline for the entire single end read dataset (SRX026839 and SRX026838) was monitored. Results for SRR094775 are given in Figure 3 since this data consists of highest number of filtered reads. Results for other datasets have been provided in Additional file 6: Figure S3. The CPU utilization in the find_juncs increases from a meagre ~12% to approximately 50% with our modified pipeline i.e. by parallelizing the find_juncs (comparing Figure 1 and Figure 3). We also observed that parallelization of find_juncs is limited to a CPU utilization of ~50% in an 8 core CPU standalone system. Increasing the extent of parallelization to more than 4 processes with such dataset in an 8 core machine, increases the time of execution thereby putting the other processes in background to halt. However, the degree of parallelization can be extended provided find_juncs is executed in a standalone system having greater number of cores. Overall, the performance gain using PVT will vary depending on the volume of the input data, availability of number of CPU cores and memory of the standalone system.

PVT outperforms TopHat

The execution time for each step involved in TopHat and PVT pipelines, run on the entire dataset was noted. As SRR094770 and SRR094775 (each chosen from the test and control data set) has larger number of filtered reads, the results for this data has been shown only (Additional file 7: Figure S4). Overall there was a significant reduction of ~23% in execution time using PVT, when compared to TopHat execution. The improvement of PVT over TopHat in find_juncs step for the entire single end read dataset is shown in Additional file 6: Figure S3.

We compared the performance (with respect to total execution time) of PVT with that of TopHat throughout the entire ‘spliced alignment’ step (shown in Additional file 8: Figure S5) for the 16 input single end read datasets (as given in Additional file 1: Table S1). Figure 4 shows that there is a strong correlation between the number of filtered reads and the total execution time in case of both TopHat and PVT. Thus we were able to predict the execution time based on the volume of the input data (maximum prediction error was within 30% of deviation). We also assessed the prediction accuracy (prediction of execution time based on the volume of the input data) of the two pipelines by fitting a linear regression model to the curve, plotting time vs. number of filtered reads as shown in Figure 4. We obtained an R-squared value of 0.87 for TopHat whereas for PVT pipeline, R-squared value improved to 0.92. In case of TopHat there is a steep increase in the execution time when the input data size exceeds the available RAM size. This may be due to the fact that TopHat does not clear garbage memory after every step of execution thereby increasing memory thrashing. However, the steepness of the slope decreases for PVT (observed with increase in R-squared value) when we modified the execution by clearing the garbage memory at every step of the pipeline. We also observed that for smaller datasets (i.e. datasets with input size less than the available RAM size) the performance gain of PVT over TopHat was low but as the number of filtered reads increases, the performance gain increases.

Using PVT for ‘spliced alignment’, the average CPU utilization increased from ~39% (as obtained from TopHat) to as high as ~66% for the entire single end read dataset (the average cache memory and CPU utilization for all the datasets used is given in Additional file 9: Figure S6 (A) and Additional file 9: Figure S6 (B) respectively).

TopHat and its pitfalls

Paired end reads (or mate pair) consist of sequence reads that occur in pairs (a constant distance of a few kb is maintained between the reads) and are obtained on sequencing both ends of the RNA simultaneously. This technique is used mainly to resolve large structural arrangements like insertions, deletions and inversions or repetitive regions of the transcriptome and is also useful for de novo transcriptome assembly. Here we modified PVT for efficient ‘spliced alignment’ of paired end reads. The dataset SRR1027730 (obtained from NCBI) for paired end reads (Additional file 2: Table S2) comprises of 45 M reads (9.5 Gbases).

The workflow of TopHat for paired end read analysis is given in Additional file 4: Figure S1. Additional file 10: Figure S7 gives the utilization graphs (CPU utilization (%) and cache memory usage) during TopHat execution of the paired end read dataset, where the steps denoted by ‘L’ are the ones for the left kept reads and ones denoted by ‘R’ are the right kept reads. The execution time for the different steps is shown in Figure 5 for SRR1027730. The total execution time for ‘spliced alignment’ of paired end reads using TopHat takes 4 hrs 44 mins for this dataset in the specified system configurations (as mentioned in Methods Section).

On analysing the execution log file for TopHat, we observed that during the execution of the steps gene_align, genome_align , junc_align and span_reads, the left kept reads and right kept reads are processed independently of each other. These independent steps can be executed simultaneously which can reduce the entire job execution time. But TopHat performs sequential execution of the independent steps along with the dependent steps, thereby increasing the waiting time for the dependent step to begin execution. This increases the overall execution time for analysing paired end read data set using TopHat.

Moreover, TopHat breaks the read into segments. Since these segments are mapped against the reference genome (during genome_align execution) and against the spliced junctions (during junc_align execution) independently, these steps can be parallelized for both the left and right kept reads. Execution time of TopHat can be further reduced by building the bowtie reference indices only once for the analysis of the entire data set of an experiment. These motivated us to design PVT for paired-end reads where we have addressed these disadvantages.

PVT

In PVT for paired end reads, we scheduled the pipeline so as to run the independent jobs simultaneously on the ‘host machine’ as well as on the ‘remote machine’ so as to reduce the entire job execution time. The workflow is shown in Additional file 11: Figure S8.

The execution of the filter_reads step generates two output files corresponding to the left and right kept reads respectively. While running TopHat, we observed that the output containing lesser number of filtered reads (either of left kept or right kept) after executing filter_reads takes lesser execution time for subsequent steps as that compared to an output having higher number of reads. Again, there was a waiting time for the dependent processes to begin their execution while running TopHat. Hence, for an efficient analysis it is evident to transfer and further process the output of the step filter_reads containing lower number of reads to the remote machine since its transfer time (both uplink and downlink) and execution time would be less than the transfer and execution of its counterpart with higher number of reads. Thereby, the execution of the steps (gene_align and genome_align; junc_align and span_reads) and required transfers for the low read files, would complete simultaneously with the execution of aforementioned steps for higher read files. This will effectively reduce the waiting time for executing the dependent steps.

We denoted the output containing lower number of filtered reads by ${\mathsf{R}}_{\mathsf{Low}}$ (coming from either left kept or right kept) and the other having higher number of filtered reads by ${\mathsf{R}}_{\mathsf{High}}$. Hence for efficient execution of subsequent steps, ${\mathsf{R}}_{\mathsf{Low}}$ is transferred to the remote machine while ${\mathsf{R}}_{\mathsf{High}}$ is executed on the host machine. After filter_reads, the execution for ${\mathsf{R}}_{\mathsf{Low}}$ and ${\mathsf{R}}_{\mathsf{High}}$ can be forked into two simultaneous alignment processes executing the consecutive steps, gene_align and genome_align.

During find_juncs, the output of genome_align from both ${\mathsf{R}}_{\mathsf{Low}}$ and ${\mathsf{R}}_{\mathsf{High}}$ are required as input to find the spliced junctions. This requires transfer of the output generated in the remote machine (on the execution of gene_align, genome_align on ${\mathsf{R}}_{\mathsf{Low}}$) to the host machine. Thus find_juncs is a dependent step, which waits for the two simultaneous processing of the outputs (i.e. ${\mathsf{R}}_{\mathsf{Low}}$ and ${\mathsf{R}}_{\mathsf{High}}$) of genome_align to complete. After the spliced junctions are found, the output for R _Low is again transferred to the remote machine for consecutive execution of the steps junc_align and span_reads. Simultaneously these steps for ${\mathsf{R}}_{\mathsf{High}}$ is executed on the host machine. In the final step report, the aligned output obtained from ${\mathsf{R}}_{\mathsf{Low}}$ and ${\mathsf{R}}_{\mathsf{High}}$ after execution of the previous steps (as mentioned above) can be concatenated, thereby generating the final alignment output and possible spliced junctions.

Using the above equations (1) and (2), we theoretically estimated the total PVT execution time for SRR1027730 based on TopHat execution time needed for each step. We obtained all the execution times $\left({\mathsf{t}}_{\mathsf{E}{\mathsf{R}}_{\mathsf{Low}}},{\mathsf{t}}_{\mathsf{E}{\mathsf{R}}_{\mathsf{High}}},{\mathsf{T}}_{\mathsf{FR}},{\mathsf{T}}_{\mathsf{FJ}}\mathsf{and}{\mathsf{T}}_{\mathsf{R}}\right)$ from TopHat execution times noted for paired end read analysis. The transfer times in PVT $\left({\mathsf{t}}_{\mathsf{T}\to {\mathsf{R}}_{\mathsf{Low}}}\mathsf{and}{\mathsf{t}}_{\mathsf{T}\leftarrow {\mathsf{R}}_{\mathsf{Low}}}\right)$ for the required inputs and outputs (for the different sets of independent steps) using scp were noted for the corresponding inputs and outputs in TopHat.

In the experiment, we implemented the PVT pipeline for SRR1027730 in two similar standalone configurations and obtained the time duration for each of the steps as given in Table 1. The improvement of PVT over TopHat as observed for the same data (shown in Additional file 12: Figure S9) experimentally was ~41% which is more than that as obtained by theoretical calculation i.e. ~34%. Such added improvement might be due to additional modification of PVT over TopHat which has been done by parallelizing sub-steps of genome_align and junc_align. This time reduction couldn’t be taken into account while calculating performance improvement of PVT over TopHat theoretically.