A GPU-accelerated metagenomic species abundance determination method

Abstract

The present invention discloses a method for determining the abundance of metagenomic species based on GPU acceleration, first analyzing the needs of user input, and determining the mode of calculation according to the input parameters of the user; then initializing, checking the integrity and version of the file, to ensure the stable operation of the present invention; after ensuring the stability, data preparation is performed, and preparation is made for subsequent GPU calculations, and all data are stored in an array of numpy format; after the current data preparation is completed, the GPU will be called by the numba package of python, and the data will be first transferred to the GPU, and then the GPU call calculation will be performed twice, the former as an intermediate result, and the latter outputs the relative abundance of each final clade; finally, the present invention will select the file format of the output according to the parameters output by the user, and save it in the saved file input by the user. The present invention greatly improves the overall efficiency of metagenomic analysis, and the accuracy is also well performed in the same type of metagenomic analysis software, with strong practical significance.

Description

Technical Field

[0001] The present invention belongs to the field of computer technology, and in particular relates to a metagenomic species abundance determination method. Background Art

[0002] A graphics processing unit (GPU), also known as a display core, graphics card, visual processor, display chip or drawing chip, is a microprocessor that specializes in performing graphics operations on personal computers, workstations, game consoles and some mobile devices. GPUs typically have hundreds or thousands of cores and are optimized to run a large number of calculations in parallel. Although GPUs are best known for 3D rendering in games, they are particularly useful for running analytics, deep learning and machine learning algorithms, allowing certain calculations to be 10 to 100 times faster than running the same calculations on traditional CPUs.

[0003] Three core ideas of GPU:

[0004] 1. Circuit simplification: The operations processed on the GPU are generally like a streaming process, without many branch conditions and complex dependencies. Therefore, the relevant circuits can be directly removed, and only the "instruction fetch", "instruction decoding", "ALU" and the registers and caches required to execute the calculations are retained.

[0005] 2. Multi-core parallelism and SIMT: Since GPU operations are naturally parallel, there is no need to use multi-threading like the CPU to achieve parallel computing. It is only necessary to integrate many groups of parallel GPU circuits in the same GPU.

[0006] 3. Hyper-threading in GPU: After GPU develops into a general-purpose computing architecture, it needs to support branch operations, but the circuits related to branch prediction have been streamlined, so it will also encounter problems similar to "pipeline stalls". GPU solves this problem through technologies similar to "hyper-threading", which is to provide more execution contexts for different tasks. Therefore, the number of execution contexts in a core is more than ALU.

[0007] GPUs optimized in the above three aspects are more adept at brute force computing. As a general computing framework, they can support different algorithms through programming. For tasks that do not have so many complex logics and conditional branches, GPUs are very suitable for computing.

[0008] Metagenome, also known as community genome, refers to the sum of all microbial genetic material in a specific environment. Metagenomics refers to the direct application of genomics technology to study microbial communities in the environment, and does not require the separation and cultivation of a single strain. At present, with the development of the times, the amount of data in metagenomics is getting larger and larger, with tens of hundreds of G at a time, resulting in more and more complex species determination of the metagenome, and longer and longer calculation time, which seriously affects the efficiency of the bioinformatics field. Moreover, the existing metagenome sequencing data analysis method is difficult to significantly speed up the analysis speed and shorten the analysis time on the basis of ensuring the accuracy of the identification results. Therefore, the bioinformatics field is in urgent need of a method to improve the computational efficiency and not affect the accuracy of the analysis results. Since the GPU can decompose a large number of operations into multiple small operations, the computational efficiency is improved by the parallel operation of a large number of threads, so the present invention uses a GPU parallelization method to accelerate the relative abundance analysis of the species of the metagenome. Summary of the invention

[0009] In order to overcome the deficiencies in the prior art, the present invention provides a method for determining the abundance of metagenomic species based on GPU acceleration, first analyzing the needs of user input, and determining the mode of calculation according to the input parameters of the user; then initializing, checking the integrity and version of the file, to ensure the stable operation of the present invention; after ensuring stability, data preparation is performed, and subsequent GPU calculations are prepared, and all data are stored in an array in numpy format; after the current data preparation is completed, the GPU will be called by the numba package of python, and the data will be first transferred to the GPU, and then the GPU call calculation will be performed twice, the former as an intermediate result, and the latter outputs the relative abundance of each final clade; finally, the present invention will select the file format of the output according to the parameters output by the user, and save it in the saved file input by the user. The present invention greatly improves the overall efficiency of metagenomic analysis, and accuracy is also well performed in the same type of metagenomic analysis software, with strong practical significance.

[0010] The technical solution adopted by the present invention to solve the technical problem includes the following steps:

[0011] Step 1: Initialization;

[0012] Read user input; check whether the bowtie2 database and metaphlan database are installed. If they are not installed or the database version is not the latest, install them automatically; and check whether bowtie2 is executable. If bowtie2 is not installed, remind the user to install it.

[0013] Read the file type input by the user. If input_type is fasta or fastq, call bowtie2 to process it and convert it to bowtie2out type. Before calling bowtie2, check the file integrity and data integrity of bowtie2.

[0014] Step 2: Data preparation;

[0015] Step 2-1: Read pkl data;

[0016] Read the pkl file in the metaphlan database. The pkl file contains the correspondence between markers and clade genes in the metagenome and stores it after reading;

[0017] Step 2: Initialize npy file;

[0018] When running the metagenomic species abundance determination process for the first time, the npy file is initialized and generated, specifically:

[0019] Convert all the data in the metagenomic species abundance determination process into numpy format, and then save the numpy format data list as an npy file;

[0020] Step 3: Preliminary data preparation;

[0021] All data calculated on the GPU are stored in numpy arrays, including: creating a correspondence between markers and reads. Each marker can correspond to multiple reads, but each read will only correspond to a unique marker; creating a parent-child relationship for each clade. Clades are divided into seven categories: kingdom, phylum, class, order, family, genus, and species. In this way, the parent-child relationship is established, and numpy arrays are used to save the parent node and all the child nodes corresponding to the parent node; finally, prepare the content required for the output file;

[0022] Step 3: GPU accelerated calculation;

[0023] Transfer the numpy array in step 2 to the GPU, and then call the GPU for calculation through the numba package; GPU calculation is divided into two parts. The first part allocates the GPU first, then calculates the removed marker corresponding to each clade, and then calculates rat_nreads based on removed, and sorts them according to the size of nreads, and finally calculates rat and loc_ab as data input for the second part; the second part will reallocate the GPU and calculate the abundance of each clade based on the output data of the previous part; finally, the data is transferred back to the CPU for integration on the CPU;

[0024] Step 4: Output;

[0025] The final output file format is determined based on the user's input parameters; before output, the metagenome is analyzed and calculated, and the final output file will be stored in the path entered by the user.

[0026] Preferably, the final output file format is CAMI or BOIM format.

[0027] Preferably, there are two methods for analyzing and calculating the metagenome in step 4: the first is rel_ab, which analyzes the metagenome from the perspective of relative abundance; the second is rel_ab_w_read_stats, which analyzes a metagenome based on relative abundance and estimates the number of reads from each branch.

[0028] The beneficial effects of the present invention are as follows:

[0029] The present invention can greatly improve the efficiency of metagenomic population analysis and abundance determination by calling GPU for multi-threaded computing, which is 2-3 times more efficient than traditional metagenomic abundance calculation software, and the efficiency of abundance calculation is improved by 6-7 times, which is a considerable breakthrough and greatly improves the overall efficiency of metagenomic analysis. In addition, after testing, the accuracy of the present invention also performs well in the same type of metagenomic analysis software, which has strong practical significance. BRIEF DESCRIPTION OF THE DRAWINGS

[0030] Figure 1 The figure is an overall flow chart of the method of the present invention. DETAILED DESCRIPTION

[0031] The present invention is further described below in conjunction with the accompanying drawings and embodiments.

[0032] The present invention provides a method for qualitative and quantitative analysis of flora based on GPU acceleration for metagenomic sequencing data. The whole process can include four continuously running modules, namely, an initialization module, a data preparation module, a GPU accelerated calculation module, and an output module. The present invention improves the efficiency of metagenomic species determination by GPU acceleration, and the calculation result is more accurate than other metagenomic species determination software. For the field of species determination in metagenomics, a new idea and method is provided.

[0033] Step 1: Initialization;

[0034] First, the user's input needs to be read. The present invention provides a large number of input parameters for the user, which can be queried through --help. After reading the user's input, for stability considerations, the present invention will also check whether the database is installed, including the bowtie2 and metaphlan databases. If they are not installed or the database version is not the latest, they will be automatically installed.

[0035] Then, the present invention will read the file type input by the user. If input_type is fasta or fastq type, bowtie2 will be called for processing and converted into bowtie2out type for subsequent processing of the program. Before calling bowtie2, the present invention will also check the file integrity and data integrity of bowtie2. Only when all the above conditions are met will the call be made to ensure the normal execution of the program.

[0036] Step 2: Data preparation;

[0037] There are three main steps:

[0038] Step 2-1: Read pkl data;

[0039] The present invention will read the pkl file in the metaphlan database, which contains the correspondence between marker and clade genes in the metagenome, and store it after reading for subsequent use.

[0040] Step 2-2: Initialize npy file;

[0041] Since the present invention needs to be calculated on the GPU, and the GPU is more friendly to the numpy data structure, all the data needed will be converted to the numpy format. In order to increase efficiency, the present invention will save these numpy format lists as npy files for use in the next execution. Therefore, when the present invention is run for the first time, it is necessary to initialize and generate npy files.

[0042] Step 2-3: Preliminary data preparation;

[0043] In order to ensure the smooth execution of the present invention on the GPU, all data that need to be calculated on the GPU must be stored in a numpy array, including: creating a correspondence between markers and reads, each marker can correspond to multiple reads, but each read will only correspond to a unique marker; creating a parent-child relationship corresponding to each clade, clades are divided into seven categories: kingdom, phylum, class, order, family, genus, and species, to establish a parent-child relationship, and using numpy to save the parent node and all the child nodes corresponding to the parent node. Finally, prepare the content required for the output file.

[0044] Step 3: GPU accelerated computing

[0045] The numpy array files prepared in step 2 are still stored on the CPU. The present invention needs to be calculated on the GPU, so these numpy array files must be transferred to the GPU, and then the GPU is called by numba for calculation. GPU calculation is divided into two parts. The first part first allocates the GPU, and then calculates the removed marker corresponding to each clade, and then calculates rat_nreads according to removed, and sorts according to the size of nreads, and finally selects the calculation method according to stat. There are 8 methods in total, and finally calculates rat and loc_ab as the data input of the second part. The second part will reallocate the GPU and calculate the abundance of each clade according to the output data of the previous part. Finally, the data is transferred back to the CPU and integrated on the CPU.

[0046] Step 4: Output;

[0047] The present invention determines the final output file format according to the user's input parameters, and provides outputs in multiple formats such as CAMI and BOIM, which are widely used in the field of bioinformatics. Before output, a step of calculation is required to determine the content of the output file. The present invention provides two commonly used calculation methods in the field of metagenomic analysis, and users can choose by themselves, and the first calculation method is used by default. The final output file will be stored in the path entered by the user. In the challenge of participating in CAMI, the present invention has F1 score = 1, False positives = 0, False negatives = 0 in the test of high-complexity data of CAMI1, which is consistent with the results of the gold standard and achieves excellent results. Specific embodiment:

[0049] The purpose of the present invention is to improve the efficiency of relative abundance calculation of metagenomes, and to use multi-threaded computing of GPU to decompose the task into multiple small tasks, and each thread performs one or more small tasks. Since the number of threads of GPU is large, and the amount of calculation of each thread is not large, the time consumption of single thread is extremely short, so as to achieve the acceleration effect. The functions of the present invention are realized by four main modules in general.

[0050] This embodiment uses Python programming language as the programming language. Figure 1In the initialization module, first, it is necessary to read various parameters input by the user, including file location parameters, request parameters, comparison parameters, other analysis type parameters and output parameters, and then save these parameters in the pars list. Then the present embodiment will check whether the metaphlan database is installed, and if it is not installed, it will automatically download and install. Next, it is necessary to set the path of the pkl data packet and bowtie2 database according to the parameters input by the user, and then check whether the pkl data packet exists. If it does not exist, an error will be reported and the program will be stopped. The file type of the input file is determined by the parameters input by the user. If the file is an unprocessed fastq or fasta type, bowtie2 will be called for processing. Before Bowtie2 processing, it will also check whether bowtie2 exists and the integrity of the bt2 file. After everything is normal, bowtie2 will be run again for processing, and a bowtie2out file will be generated and saved. If the same fastq file is input next time, the present embodiment will remind the user that the fastq file has generated a bowtie2out file, and the file can be directly used as an input file to speed up the operation efficiency.

[0051] Figure 1In the data preparation module, it is necessary to read bowtie2out or sam file first, which includes the mapping relationship between marker and reads. Then read the content in the pkl data packet, which contains all marker information in metaphlan. After the file is read, it will be checked whether this embodiment is running for the first time. If it is running for the first time, it is necessary to prepare npy files for GPU calculation. The marker id, name information, the correspondence between marker and extid, the correspondence between marker and clade, the father-child correspondence of clade, the id and name information of clade, and the array information of intermediate variables such as clade_t and vir_cladeid are stored by the savenpy() function. The above data will only be generated when the present embodiment is run for the first time. If it has been generated before, it will not be repeated, so as to speed up the efficiency of the execution of the present embodiment. Read the above npy file for subsequent GPU calculation. Read the bowtie2out file into the dictionary marker2reads, and then use the above data to generate the correspondence between clade and marker. Each clade will correspond to a marker, and each marker can correspond to multiple clades. In addition, each clade marker will correspond to multiple reads, so the TaxClade_markers2nreads dictionary is needed to store the number of reads corresponding to each clade marker. All data is stored in numpy format for subsequent GPU use.

[0052] Figure 1In the GPU accelerated calculation module, this embodiment uses the numba package in python to call the GPU, which is more in line with the programming habits of the python language than cuda programming. In this module, due to the many restrictions on GPU calculations, part of the calculations cannot be run on the GPU, so it needs to be divided into two parts and put on the GPU for calculation, and some subsequent operations are performed on the CPU. Since the storage of the GPU and the CPU is not common, it is necessary to transfer the data from the CPU to the GPU through the cude.to_device() function before calculating, and then the threads can be allocated, and each thread performs one or more simple tasks. In the first GPU calculation, the nreads information of the marker corresponding to each clade is stored. This embodiment also provides 8 clade processing methods, and users can choose according to their needs. After the GPU calculation, two intermediate data arrays of locab and rat will be output. Next, the program will execute cuda.synchronize() to wait for all threads in the GPU to complete execution, and then proceed to the next step of GPU calculation. According to the seven methods of kingdom, phylum, class, order, family, genus and species, they are imported into the GPU for calculation, and the GPU is used to accelerate the summation to obtain the abundance information of each clade. Finally, the data needed for CPU calculation is copied from the GPU to the CPU through the copy_to_host() function, and then sorting and summing operations are performed.

[0053] Figure 1 In the output module, this embodiment needs to first determine the header content of the output file according to the user's input parameters, and provides a variety of output formats, the default output format of metaphlan. After writing the header content, it is also necessary to analyze the results of the GPU accelerated calculation module. This embodiment provides two analysis and calculation methods. The first is "rel_ab", which is also the default calculation and analysis method. This method will analyze the metagenome from the perspective of relative abundance. The second is "rel_ab_w_read_stats", which analyzes a metagenome based on relative abundance and estimates the number of reads from each branch. These two methods are commonly used analysis methods in metagenomic analysis. The two methods are the same in GPU distributed computing, but there are some differences in the early data preparation and subsequent output content.

Claims

1. A GPU-accelerated metagenomic species abundance determination method, characterized in that: The steps include: Step 1: Initialization; Read user input; check whether the bowtie2 database and metaphlan database are installed. If they are not installed or the database version is not the latest, install them automatically; and check whether bowtie2 is executable. If bowtie2 is not installed, remind the user to install it. Read the file type input by the user. If input_type is fasta or fastq, call bowtie2 to process it and convert it to bowtie2out type. Before calling bowtie2, check the file integrity and data integrity of bowtie2. Step 2: Data preparation; Step 2-1: Read pkl data; Read the pkl file in the metaphlan database. The pkl file contains the correspondence between markers and clade genes in the metagenome and stores it after reading; Step 2: Initialize npy file; When running the metagenomic species abundance determination process for the first time, the npy file is initialized and generated, specifically: Convert all the data in the metagenomic species abundance determination process into numpy format, and then save the numpy format data list as an npy file; Step 3: Preliminary data preparation; All data calculated on the GPU are stored in numpy arrays, including: creating a correspondence between markers and reads. Each marker can correspond to multiple reads, but each read will only correspond to a unique marker; creating a parent-child relationship for each clade. Clades are divided into seven categories: kingdom, phylum, class, order, family, genus, and species. In this way, the parent-child relationship is established, and numpy arrays are used to save the parent node and all the child nodes corresponding to the parent node; finally, prepare the content required for the output file; Step 3: GPU accelerated calculation; Transfer the numpy array in step 2 to the GPU, and then call the GPU for calculation through the numba package; GPU calculation is divided into two parts. The first part allocates the GPU first, then calculates the removed marker corresponding to each clade, and then calculates rat_nreads based on removed, and sorts them according to the size of nreads, and finally calculates rat and loc_ab as data input for the second part; the second part will reallocate the GPU and calculate the abundance of each clade based on the output data of the previous part; finally, the data is transferred back to the CPU for integration on the CPU; Step 4: Output; The final output file format is determined based on the user's input parameters; before output, the metagenome is analyzed and calculated, and the final output file will be stored in the path entered by the user.

2. The GPU-accelerated metagenomic species abundance determination method according to claim 1, characterized in that: The final output file format is CAMI or BOIM format.

3. The method for determining species abundance based on GPU acceleration based on metagenomics according to claim 1, characterized in that: There are two methods for analyzing and calculating the metagenome in step 4: the first is rel_ab, which analyzes the metagenome from the perspective of relative abundance; the second is rel_ab_w_read_stats, which analyzes a metagenome based on relative abundance and estimates the number of reads from each branch.

Images