JP2020509473A - Compact representation method and apparatus for biological information data using a plurality of genome descriptors - Google Patents

Abstract

A method and apparatus for compressing genomic sequence data generated by a genomic sequencing device. Sequence reads are coded by aligning against an existing or constructed reference sequence, the process of coding being based on the classification of reads into data classes, followed by multiplexed descriptor blocks. It consists of coding for each class. A specific source model and entropy coder is used for each data class into which the data is divided, and each associated descriptor block. [Selection diagram] FIG.

Description

The present disclosure provides a new representation of genomic sequencing data that provides new functionality not available with known prior art representations, thereby reducing storage space used and improving access performance.
[Cross-reference of related applications]

  This application claims the benefit of PCT / US2017 / 017842, filed February 14, 2017 and PCT / US2017 / 041591, filed July 11, 2017, and their benefits.

  Proper representation of genomic sequencing data is needed to enable efficient genomic analysis applications such as genomic variant calls and all other analyzes performed for various purposes by processing sequence data and metadata. It is essential.

  Sequencing of the human genome has become cheaper with the advent of high-throughput, low-cost sequencing technology. Such opportunities range from pathogen surveillance for the identification of antibodies, from cancer diagnosis and treatment to the identification of genetic disorders, to the creation of new vaccines, drugs, and the customization of personalized treatments. It opens up new perspectives in several areas.

  Hospitals, genomics data analytics providers, bioinformatics, and large biological data storage centers are inexpensive, fast, reliable, and interconnected that enable genomic medicine to scale up globally Looking for a genome information processing solution. One of the bottlenecks in the sequencing process is data storage, and methods of expressing genomic sequencing data in a compressed format are being increasingly studied.

  The display of genomic information most used in sequencing data is based on FASQ and SAM format compression. The purpose is to compress the conventionally used file formats (for unaligned data and aligned data use FASTQ and SAM respectively). Such files are composed of plain text characters and described above using a generic approach such as the LZ (Lempel and Ziv, author of the first edition) (well-known zip, gzip, etc.) Compressed. When using a general purpose compression scheme such as gzip, the result of the compression is typically a single chunk of binary data. Such monolithic information becomes extremely difficult to archive, transfer, and refine, especially when the amount of data is very large, as in the case of high-throughput sequencing. The BAM format is inefficient, rather than extracting the actual genomic information conveyed by the SAM file, and focuses on the compression of the redundant SAM format, and rather than exploiting the specific properties of each data source It is characterized by low compression performance in order to adopt a general purpose text compression algorithm such as gzip (genomic data itself).

  A less sophisticated but more sophisticated approach to genomic data compression that is more efficient than BAM is CRAM. CRAM provides efficient compression by employing differential coding on the reference (partially exploits data source redundancy). However, features such as incremental updates, streaming support, and selective access to certain classes of compressed data are not yet available.

  With these approaches, the compression ratio is low and the data structures are compressed, which makes navigation and manipulation difficult. Even when performing simple operations or accessing selected regions of a genomic dataset, downstream analysis can be very slow due to the need to process large and rigorous data structures . CRAM relies on the concept of CRAM records. Each CRAM record represents a single mapped or unmapped lead by encoding all the elements required for reconstruction.

CRAM has the following disadvantages and limitations that are solved and overcome by the invention described herein:
1. CRAM does not support random access to data indexes and data subsets that share certain functions. The index of the data is outside the scope of the specification (see section 12 of CRAM specification v.3.0) and implemented as a separate file. In contrast, the inventive approach described herein employs a data indexing method that is integrated with the encoding process, with the index embedded in the encoded (ie, compressed) bitstream. I have.
2. A CRAM is constructed with a core data block that can contain any type of mapped leads (exact matches, leads with only replacements, leads with insertions or deletions (also called "indels")). Is done. There is no concept of classifying data or grouping leads within classes according to the result of mapping for the reference sequence. This means that all data must be inspected, even if only leads with a particular function are searched. Such limitations are solved by the present invention, which classifies and divides data by class before coding.
3. CRAM is based on the concept of encapsulating each lead into a "CRAM record". This means that when searching for reads characterized by a particular biological feature (eg: reads with substitutions but without “indels” or fully mapped reads), the full version of each will be searched. It means the need to inspect the "record".
In contrast, in the present invention, there is the concept of a separately coded data class in a separate block of information and no record concept to encapsulate each lead. This allows for the reading of specific biological properties (eg: reads with substitutions but without "indels", or fully mapped More efficient access to the set of leads having the same read).
4. In a CRAM record, each record field is associated with a particular flag, and since each CRAM record can contain different types of data, there is no concept of a context and each flag must always have the same meaning. This coding mechanism introduces redundant information and precludes the use of efficient context-based entropy coding.
On the other hand, in the present invention, the flag indicating the data is essentially defined by the "block" of the information to which the data belongs, so that there is no concept of the flag indicating the data. This means that the number of symbols to be used is significantly reduced, and consequently the entropy of the information source, which results in more efficient compression. Such an improvement is possible because the use of different "blocks" allows the encoder to reuse the same symbols across blocks that have different meanings depending on the context. In a CRAM, since there is no concept of a context and each CRAM record can include any kind of data, each flag must always have the same meaning.
5. In CRAM replacement, insertions and deletions are represented using different descriptors, options that increase the size of the alphabet of the information source and result in higher entropy of the information source. In contrast, the disclosed inventive approach uses a single alphabet and coding for substitutions, insertions and deletions. This simplifies the coding and decoding process, yields a bitstream in which the coding is characterized by high compression performance, and produces a source model with low entropy.

  The present invention classifies the sequencing data so that redundant information to be encoded is minimized and features such as support for selective access and incremental updates are directly possible within the compression domain. And compressing the genome sequence by dividing.

  The following features according to the claims solve the problems of existing prior art solutions by their provision.

A method for encoding genomic sequence data comprising a nucleotide sequence read, comprising:
Aligning said lead with one or more reference sequences, thereby creating an alignment lead;
Classifying said alignment reads using said one or more reference sequences according to a specified matching rule, thereby creating a class of alignment reads;
Encoding the classified alignment reads as multiple blocks of descriptors,
Encoding the classified alignment reads as multiple blocks of the descriptor includes selecting the descriptor according to the class of the alignment reads,
The block of descriptors is structured with header information, thereby creating a continuous access unit.

In another aspect, the coding method further comprises classifying the leads that do not satisfy the specified matching rule into a class of unmapped leads,
Constructing a set of reference sequences using at least some of the unmapped reads;
Aligning the unmapped class of reads with the constructed set of reference sequences,
Encoding the classified alignment reads as multiple blocks of descriptors,
Encoding the set of constructed reference sequences;
The block of descriptors and the coded reference sequence are constructed with header information, thereby creating successive access units.

  In another variation, the encoding method further comprises classifying genomic reads without mismatches in the reference sequence as a first “Class P”.

  In another embodiment, the encoding method finds a mismatch only at locations where the sequencing device cannot call any "bases" and the number of mismatches in each read does not exceed a predetermined threshold. In some cases, further comprising classifying the genomic read as a second "Class N".

  In another embodiment, the encoding method comprises identifying a genomic read as a third "Class M" if a mismatch is found at a position where the sequencing device was unable to call any "bases". Further referred to as "n-type" mismatches and / or referred to as "bases" different from the reference sequence, referred to as "s-type" mismatches, and wherein the number of mismatches is greater than the "n-type" mismatches. Does not exceed a predetermined threshold for the number of "s-type" mismatches, and the threshold is a function (f (f ()) that calculates the number of "n-type" and "s-type" mismatches. n, s)).

  In another embodiment, the encoding method further comprises identifying the genomic read as a fourth "Class I" if a mismatch of the same type as the "Class M" is likely to occur, at least One type of mismatch: “insert” (“i-type”), “delete” (“d-type”), soft clip (“c-type”), where the number of mismatches of each type is Without exceeding the corresponding predetermined threshold, the threshold is a function (w) that calculates the number of "n-type", "s-type", "i-type", "d-type" and "c-type" mismatches. (N, s, i, d, c)).

  In another embodiment, the encoding method comprises identifying the genomic read as a fifth "class U" as including all reads that do not find any of the classes P, N, M, I. In addition.

  In another embodiment, the encoding method further comprises that the encoded reads of the genomic sequence are paired.

  In another aspect, the method of encoding, wherein the classifying comprises all read pairs where one read belongs to class P, N, M or I and the other read belongs to "class U", The method further includes identifying the genomic read as a sixth "Class HM".

In another aspect, the encoding method identifies whether the reads of the two mates are classified into the same class (P, N, M, I, U, respectively), and identifies the pair as the same. Assigned to class,
Identify whether the two mate's leads are classified into different classes, and if neither belongs to "Class U", assign the paired leads to the highest priority class according to:
P <N <M <I
Here, "class P" has the lowest priority, "class I" has the highest priority,
Identifying whether only one of the two mates' leads is classified as belonging to "Class U" and classifying the paired leads as belonging to a "Class HM" sequence.

  In another aspect, the encoding method further comprises: wherein each class of leads N, M, and I has a number (292) of “n-type” mismatches, a function f (n, s) (293), and a function w (n, s, i, d, c) (294) and two or more subclasses (296, 297) according to the threshold vectors (292, 293, 294) defined for each class N, M, I, respectively. , 298).

Identifying whether the leads of the two mates are classified in the same subclass, assigning the pair to the same subclass,
Identifying whether the leads of the two mates are classified into different classes of subclasses and assigning the pair to the subclass belonging to the higher priority class according to the following formula:
N <M <I
Here, N has the lowest priority, I has the highest priority,
Identify whether the two mates' leads are classified into the same class, which class is N, M, or I, but have different subclasses, and assign the pair to the highest priority according to the following formula: Assign to higher subclass,
N ₁ <N ₂ <... <N _k
M ₁ <M ₂ <... M _j
I ₁ <I ₂ <... <I _h
Here, it further includes that the highest index has the highest priority.

  In another aspect, information regarding the mapping location of each lead is encoded by a pos descriptor block.

  In another aspect, information about the strandiness of each read (ie, the DNA strand from which the sequence of reads is derived) is encoded by an rcomp descriptor block.

  In another aspect, the pairing information of the paired-end read is encoded by a pair descriptor block.

  In another embodiment, additional readings such as whether the reads are mapped in proper pairs, failure of platform / vendor quality check, PCR or optical replication, or auxiliary alignment. Alignment information is encoded by a flag descriptor block.

  In another aspect, the information about the unknown base is encoded by a mmis descriptor block.

  In another aspect, the information regarding the location of the substitution is encoded by a snpp descriptor block.

  In another aspect, the information regarding the type of substitution is encoded by a particular snpt descriptor block.

  In another aspect, the information about the location of the mismatch, substitution, insertion or deletion is encoded by an indp descriptor block.

  In another aspect, information about the type of mismatch, such as a substitution, insertion, or deletion, is encoded by an indt descriptor block.

  In another aspect, information about the clipped bases of the mapped read is encoded by an indc descriptor block.

  In another aspect, information about unmapped leads is encoded by a uleads descriptor block.

  In another aspect, information about the type of reference sequence used for coding is coded by an rtype descriptor block.

  In another aspect, the information regarding multiple alignment of the mapped reads is encoded by an mmap descriptor block.

  In another aspect, the information regarding the spliced alignment and the multiple alignment of the same read is encoded by an msar descriptor block and an mmp descriptor block.

  In another aspect, information about the alignment score of the read is encoded by an mscore descriptor block.

  In another aspect, information about the group to which the lead belongs is coded by an “rgroup” descriptor block.

  In another aspect, the coding method further comprises the block of descriptors includes a master index table including one section for each class and subclass of aligned reads, the section comprising a master index table and the access table. It further includes the mapping location on the one or more reference sequences of the first read of each access unit of each class or subclass of data that is encoded in both units.

  In another aspect, the encoding method further comprises information about the type of reference (existing or constructed) in which the block of the descriptor is used and the segment of the lead that is not mapped to the reference sequence. .

  In another aspect, the encoding method further comprises the reference sequence is first transformed to a different reference sequence by applying substitution, insertion, deletion, and clipping, and the classification as a number of blocks of descriptors. The coding of the alignment read further includes referencing the converted reference sequence.

  In another aspect, the coding method further comprises applying the same transformation to the reference sequence used for all classes of data.

  In another aspect, the encoding method further comprises that a different transform is applied to the reference sequence used for each class of data.

  In another aspect, the coding method further comprises the transform of the reference sequence is coded as a block of descriptors and structured with header information, thereby creating a continuous access unit.

  In another aspect, the encoding method further comprises: the encoding of the classified alignment read and the associated reference sequence transformation as a multiplex of descriptor blocks comprising a particular descriptor block and a particular source model. Further comprising the step of:

  In another aspect, the coding method further includes that the entropy coder is one of a context adaptive arithmetic coder, a variable length coder, and a Golomb coder.

The present invention further provides a method of decoding encoded genomic data, comprising:
Analyzing the access unit containing the encoded genomic data to extract a multiplexed block of descriptors using the header information;
A method is provided that includes decoding the multiplexed block of descriptors to extract leads according to specific matching rules that define a classification for one or more reference sequences.

  In another aspect, the decoding method further comprises decoding the unmapped genomic reads.

  In another aspect, the decoding method further comprises decoding the classified genomic reads.

  In another aspect, the decoding method further includes decoding a master index table that includes one section for each class of associated associated mapping locations and leads.

  In another aspect, the decoding method further comprises decoding information related to the type of reference used: existing, transformed, or constructed.

  In another aspect, the decoding method further comprises decoding information related to one or more transforms applied to the existing reference sequence.

  In another embodiment, the decoding method further comprises paired genomic reads.

  In another aspect, the decoding method further includes a case where the genomic data is entropy decoded.

The present invention further provides a genome encoder (210) for compressing the genome sequence data 209 including a genome sequence data 209 and a nucleotide sequence read, wherein the genome encoder (210) comprises:
An aligner unit (201) configured to align the lead with one or more reference sequences, thereby creating an alignment lead;
A constructed reference generation unit (202) configured to generate the constructed reference sequence;
Data configured to classify said alignment reads according to a particular matching rule using one or more existing or constructed reference sequences, thereby creating a class of alignment reads (208). A classification unit (204);
One or more block coding units (205-207) configured to code the classified alignment read as a block of descriptors by selecting the descriptor according to the classified alignment read; ,
A multiplexer (2016) for multiplexing the compressed genomic data and metadata.

  In another aspect, the genomic encoder further includes a reference sequence conversion unit (2019) configured to convert an existing reference and data class (208) to a converted data class (2018).

  In another aspect, the genomic encoder further comprises an encoder of data classes N, M and I, wherein said data classification unit (204) is constituted by a vector of thresholds generating subclasses of data classes N, M and I. Including.

  In another aspect, the genomic encoder further comprises the reference transform unit (2019) applying the same reference transform (300) to all classes and subclasses of data.

  In another aspect, the genomic encoder further comprises the reference transform unit (2019) applying different reference transforms (301, 302, 303) to different classes and subclasses of data.

  In another aspect, the genomic encoder further comprises functions suitable for performing all of the above-described encoding methods.

The present invention further provides a genome decoder (218) for decompressing the compressed genome stream (211), wherein the genome decoder (218) comprises:
A demultiplexer (210) for demultiplexing the compressed genomic data and metadata;
Analysis means (212-214) configured to parse the compressed genome stream into genome blocks (215) of descriptors;
One or more block decoders (216-217) configured to decode genomic blocks of the descriptor into sorted reads of a sequence of nucleotides (211);
A genomic data class decoder configured to selectively decode said categorized reads of the sequence of nucleotides on one or more reference sequences to generate an uncompressed read of the sequence of nucleotides; ,including.

  In another aspect, the genomic decoder (2112) decodes the reference translation descriptor (2112) and generates a translated reference (2114) for use by the genomic data class decoder (219). 2113).

  In another aspect, the genomic decoder further comprises that the one or more reference sequences are stored in a compressed genomic stream (211).

  In another aspect, the genomic decoder further comprises the one or more reference sequences being provided to the decoder via an out-of-band mechanism.

  In another aspect, the genomic decoder further comprises the one or more reference sequences being constructed at a decoder.

  In another aspect, the genomic decoder further comprises that one or more reference sequences are converted at the decoder by a reference conversion decoder (2113).

  The invention further provides a computer readable medium comprising instructions for executing at least one processor for performing all aspects of the coding method described above.

  The invention further provides a computer readable medium comprising instructions for executing at least one processor for performing all aspects of the decoding method described above.

  The present invention further provides support data for storing a genome encoded according to all aspects of the encoding method described above.

One aspect of the proposed approach is the definition of separately coded data and metadata classes that are structured in different blocks. More appropriate improvements of such methods over existing methods are as follows:
1. Improving compression performance by reducing the entropy of information sources configured by providing an efficient source model for each type of data or metadata;
2. Possibility of selectively accessing compressed data and some of the metadata directly in the compressed domain for further processing purposes;
3. Include compressed data and metadata incrementally (i.e., no decoding and recoding required) with new sequencing data and / or metadata and / or new analysis results associated with a particular set of sequence reads ) Possibility to update.

FIG. 4 shows how the position of the mapped lead pair is coded in the “pos” block as the difference from the absolute position of the first mapped lead. Figure 2 shows how two reads in a pair are generated from two DNA strands. If strand 1 is used as a reference, it shows how the reverse complement of lead 2 is coded. The four possible combinations of reads that make up the read pair and their respective encodings in the "rcomp" block are shown. A method for calculating a pairing distance when the read lengths of three read pairs are constant is shown. Pairing errors coded in "pair" blocks show how the decoder reconstructs correct lead pairing using coded "MPPPD". FIG. 4 shows the encoding of the pairing distance when a lead is mapped to a different reference than its mate. In this case, an additional descriptor is added to the pairing distance. One is a signaling flag, the second is a reference identifier, and a pairing distance. FIG. 7 illustrates the encoding of “n-type” mismatches in the “nmis” block. 5 shows a mapped read pair indicating replacement for a reference sequence. A method of calculating the position of the replacement as an absolute value or a differential value will be described. The calculation method of the symbol which codes the type of substitution when the IUPAC code is not used is shown. The symbol represents the distance between the molecule present in the lead and the molecule present in the reference at that position-the circular substitution vector. Figure 3 shows how the permutation is coded into a "snpt" block. A method of calculating a replacement code when an IUPAC ambiguity code is used will be described. The coding method of the "snpt" block when the IUPAC code is used is shown. The permutation vector used in the class I lead is the same as the class M, and shows that the special codes are added to the insertion of the symbols A, C, G, T, and N. FIG. 6 shows examples of mismatch and indels coding for IUPAC ambiguity codes. In this case, the permutation vector will be very long, so there are more possible calculated symbols than with five symbols. Figure 3 shows different source models of mismatches and indels, where each block contains the location of a single type of mismatch or insertion. In this case, the symbols are not coded for mismatch or indels types. Figure 2 shows examples of mismatch and indels coding. If a particular type of mismatch or indel is not present in the lead, a 0 is coded in the corresponding block. 0 functions as a separator and a terminator for each block. The change of the reference sequence shows how to convert the M-read into the P-read. This operation can reduce the information entropy of the data structure, especially in the case of high coverage data. 4 shows a genome encoder 2010 according to one embodiment of the present invention. 4 shows a genome decoder 218 according to one embodiment of the present invention. FIG. 4 illustrates how to build an “internal” reference by clustering leads and assembling the segments obtained from each cluster. 7 illustrates how a reference is constructed by storing the latest leads after a particular sort (eg, lexicographic order) has been applied to the leads. Fig. 7 illustrates how leads belonging to the class of "unmapped" leads (class U) are encoded using six descriptors stored or transmitted in corresponding blocks. 4 shows an alternative encoding of a lead belonging to class U. Here, the mapping position of the constructed reference lead is encoded using a coded pos descriptor. 4 illustrates a method of applying a reference to remove mismatches from a read. In some cases, the reference transformation may create a new mismatch or change the type of mismatch found when referencing the reference before the transformation is applied. If all or a subset of the mismatches are deleted (ie, leads belonging to class M before conversion are assigned to class P after reference conversion is applied), how the reference conversion belongs to class lead Indicates whether the destination can be changed. FIG. 4 illustrates how to fill unknown regions of a reference sequence by constructing long contigs with unmapped reads using a half-mapped read pair (class HM). FIG. 4 shows how the encoders of the classes N, M and I data are composed of a vector of thresholds and generate separate subclasses of the N, M and I data classes. All classes of data can use the same transformed reference for recoding, or use different transformations for each class N, M and I or any combination thereof Indicate if you can. 3 shows the structure of a genome data set header. 1 shows a general structure of a master index table. Each row contains genomic intervals for several classes of data P, N, M, I, U, HM, and pointers to metadata and annotations. The columns indicate a particular location on the reference sequence associated with the encoded genomic data. FIG. 4 shows an example of a single line of an MIT containing genomic intervals associated with class P reads. Genomic regions associated with different reference sequences are separated by special flags ("S" in the example). Fig. 3 illustrates the general structure of a local index table (LIT) and the method used to store a pointer to the physical location of the encoded genomic information contained in the stored or transmitted data. 9 shows an example of an LIT used to access access unit numbers 7 and 8 of a block payload. 3 shows a functional relationship between a plurality of lines of MIT and LIT included in a genome block header. It shows how an access unit is composed of several blocks of genomic data transmitted by different genomic streams containing data belonging to different classes. Each block is further constituted by a data packet used as a data transmission unit. It shows how the access unit is composed of multiplexed blocks belonging to one or more blocks of the same data as the header. Each block can consist of one or more packets containing the actual descriptors of the genomic information. 2 shows multiple alignment without splicing. The leftmost lead has N alignments. N is the first value of mmap to be decoded, and indicates the number of alignments of the first read. The next N values in the mmap descriptor are decoded and used to calculate P, the number of alignments in the second read. FIG. 4 illustrates how to encode multiple alignments without splices using position, pair, and mmap descriptors. The leftmost lead has N alignments. 3 shows a multiple alignment using a splice. Figure 4 illustrates how pos, pair, mmap, and msar descriptors are used to represent multiple alignments with splices.

  Genomic or proteomic sequences according to the present invention include, for example, without limitation, nucleotide sequences, deoxyribonucleic acid (DNA) sequences, ribonucleic acid (RNA), and amino acid sequences. Although the description herein is fairly detailed with respect to genomic information in the form of nucleotide sequences, as will be appreciated by those skilled in the art, there are several variations, and methods and systems for compression can be performed on other genomes or It will be appreciated that the same applies to proteome sequences.

  Genomic sequencing information is generated by high-throughput sequencing (HTS) equipment in the form of sequences of nucleotides (also called "bases") represented by strings from a defined vocabulary. The minimal vocabulary is represented by five symbols: {A, C, G, T, N} represent the four nucleotides present in DNA: adenine, cytosine, guanine, thymine. Thymine is replaced by uracil (U) in RNA. N indicates that when the sequencing device failed to call any base, the actual nature of that position has not been determined. If the IUPAC ambiguity code is employed by the sequencing device, the alphabet used for the symbols is (A, C, G, T, U, W, S, M, K, R, Y, N, D, H, V, N).

  The nucleotide sequence generated by the sequencing device is called a "read". Sequence reads have a length of tens to thousands of nucleotides. In some techniques, one read generates a “pair” of sequence reads from one DNA strand and a second read from the other strand. Genome sequencing uses the term "coverage" to describe the level of redundancy in sequence data with respect to a "reference sequence." For example, to achieve 30-fold coverage in the human genome (333.2 billion bases in length), the sequencing equipment would generate a total of 30 x 3.2 billion bases and each reference position would be "covered" on average 30 times. I do.

  Throughout this disclosure, a reference sequence is any sequence to which a nucleotide sequence generated by a sequencing device is aligned / mapped. One example of a sequence is actually a “reference genome”, which is a sequence assembled by scientists as representative of a species gene set. For example, GRCh37, the human genome of the Genome Reference Consortium (build37), is derived from 13 anonymous volunteers in Buffalo, NY. However, the reference sequence can also be composed of a synthesized sequence that has been devised and constructed to simply improve in consideration of further processing of the read compressibility. This is described in detail in "Building Class U Descriptors and" Internal "References for Unmapped Leads of" Class U "and" Class HM ", and is shown in FIGS.

In the sequencing device, the following sequence read error may occur.
1. Decision to skip calling a base because there is no reliability calling a particular base. This is called the unknown base and is labeled "N" (indicated as an "n-type" mismatch).
2. Use the wrong symbol (ie, representing a different nucleic acid) to represent the nucleic acid actually present in the sequenced sample; this is commonly referred to as a "substitution error" (indicated as an "s-type" mismatch) ).
3. Insertion of one of the additional symbols into the sequence read without reference to the nucleic acid actually present; this is commonly referred to as an "insertion error" (indicated as an "i-type" mismatch).
4. Deletion of one of the symbols representing nucleic acids actually present in the sequenced sample from the sequence read; this is commonly referred to as a "deletion error" (indicated as a "d-type" mismatch).
5. Recombination of one or more fragments into a single fragment that does not reflect the reality of the original sequence; this usually results in the aligner deciding to clip the base (indicated as a "c-type" mismatch).

The term "coverage" is used in the literature to quantify the extent to which a reference genome or a part thereof can be covered by available sequence reads. Coverage is said to be:
• partial (less than 1x) where some parts of the reference genome are not mapped by any readable sequence;
A single (1x) in which all nucleotides of the reference genome are mapped by a single symbol in the sequence;
-Multiple (2x, 3x, Nx) if each nucleotide of the reference genome is mapped multiple times.

  An object of the present invention is to define a genome information display format in which related information is efficiently accessible and movable, and the weight of redundant information is reduced.

The main innovative aspects of the disclosed invention are as follows.
1 Sequence reads are classified and classified into data classes according to the result of alignment with respect to the reference sequence. Such classification and segmentation allows for selective access to the encoded data according to criteria relating to alignment results and matching accuracy.
2 The categorized sequence reads and associated metadata are represented by like blocks of descriptors to obtain a separate information source characterized by low information entropy.
3 Possibility to model each separated information source using a separate source model adapted to the statistical properties of each class, and data units (access units) within each class of the lead and each separately accessible Possibility to change the source model within each descriptor block. Adopting an appropriate context adaptive stochastic model and associated entropy coder according to the statistical properties of each source model.
4. Defining correspondences and dependencies between descriptor blocks allows for selective access to sequencing data and associated metadata without decoding all descriptor blocks when not all information is needed.
5 The data class of each sequence and the code of the associated metadata block for the “existing” (also called “external”) reference sequence or the “transformed” reference sequence are used to reduce the entropy of the information source of the descriptor block. , Obtained by applying the appropriate transformation to the "existing" reference sequence. The descriptor represents a lead divided into different data classes. Referencing the "existing" reference or the "transformed""existing" reference sequence, following the coding of the read using the corresponding descriptor, using the occurrence of various mismatches, the low entropy To find the final coded representation and achieve higher compression efficiency, an appropriate conversion to the reference sequence can be defined.
6 The construction of one or more reference sequences (also referred to as “external” reference sequences, also referred to herein as “external” reference sequences, as distinguished from “existing” reference sequences) is based on existing Used to encode a class of reads that exhibits some degree of matching accuracy with respect to the reference sequence. Such constraints imply that the cost of coding to represent the class of reads aligned with respect to the "internal" reference sequence in a compressed format, and the cost of representing the "internal" reference sequence itself, is aligned. The purpose is to set the class of unread leads to be lower than using verbatim coding or using an "external" reference sequence with or without conversion.

Hereinafter, each embodiment will be described in more detail.
[Classification of sequence read according to matching rule]

  Sequence reads generated by the sequencing device are classified into six different "classes" according to the disclosed invention, according to the alignment matching results for one or more "existing" reference sequences.

When aligning a DNA sequence of nucleotides to a reference sequence, the following cases can be identified:
1. Certain regions within the reference sequence are found to be consistent with error-free sequence reads (ie, perfect mapping). Such sequences of nucleotides are referred to as "perfectly matched reads" or "class P" Is displayed.
2. Certain regions of the reference sequence are found to be consistent with sequence reads with the number and type of mismatches determined solely by the number and location at which the sequencing device generating the read was unable to call the base (or nucleotide). . Such type of mismatch is indicated by the letter "N" used to indicate an undefined nucleotide base. This type of mismatch is referred to herein as an "n-type" mismatch. Such a sequence belongs to a "Class N" lead. When a lead is classified as belonging to "Class N", it is convenient to limit the degree of inaccuracy of the matching to a certain upper limit and to set a boundary between what is considered a valid match and what is not. Therefore, reads assigned to class N are also constrained by setting a threshold (MAXN) that defines the maximum number of undefined bases (bases called "N") that can be included in the read. . Such a classification implicitly defines the minimum matching accuracy (or maximum matching degree) that all leads belonging to class N need to share when referring to the corresponding reference sequence, which is optional. Configure a useful criterion for applying data retrieval to compressed data.
3. One region in the reference sequence is the number of positions where the sequencing device that generated the read could not call any nucleotide bases, if any (ie, an “n-type” mismatch), plus It can be seen that the number of mismatches called bases different from those present are consistent with the sequence read with the number and type of mismatches determined by the mismatch. Such a type of mismatch, denoted as "substitution", is also called a single nucleotide mutation (SNV) or a single nucleotide polymorphism (SNP). This type of mismatch is referred to herein as an "s-type" mismatch. The sequence read is referred to as “M mismatch read” and is assigned to “class M”. As in the case of "Class N", for all leads belonging to "Class M", the degree of matching inaccuracies is limited to a certain upper limit, and the difference between what is considered a valid match and what is not. It is convenient to set boundaries. Thus, the leads assigned to class M are constrained by defining a set of thresholds, one for the number “n” (MAXN) if there is an “n-type” mismatch, and one for Is the number of substitutions “s” (MAXS). A third constraint is a threshold defined by a function f (n, s) of both numerical values "n" and "s". Such a second constraint makes it possible to generate a class with an upper bound on matching inaccuracies according to any meaningful selective access criteria. For example, but not by way of limitation, f (n, s) bounds to (n + s) 1/2, or (n + s), or the maximum matching inaccuracy level allowed for leads belonging to "Class M". Can be any linear or non-linear expression that sets Such boundaries exceed the simple threshold applied to one or the other type, and further limit the possible combinations of the number of “n-type” mismatches and “s-type” mismatches (permutations). To constitute a very useful criterion for applying the desired selective data retrieval to the compressed data when analyzing sequence reads for various purposes.
4. The fourth category is composed of sequencing reads that show at least one mismatch of any type of “insert”, “delete” (also called indels), “clip”, and furthermore class N or M There is a type of mismatch belonging to Such a sequence is called "I mismatch read" and is assigned to "Class I". Insertions are made up of additional sequences of one or more nucleotides that are not present in the reference but are present in the read sequence. This type of mismatch is referred to herein as an "i-type" mismatch. When the inserted sequence is at the end of the sequence, it is also referred to in the literature as a “soft clip” (ie, the nucleotides do not match the reference, but are aligned with the discarded “hard clip” nucleotides) In the read lead). This type of mismatch is referred to herein as a "c-type" mismatch. The retention or discarding of nucleotides is determined by a sequencing device or by the aligner stage rather than by the read classifier disclosed and disclosed in the present invention that receives and processes reads, as determined by the following sequencing stages. It is. This is a decision made by the aligner stage, rather than by the sequencing device or by the lead classifier disclosed and disclosed in the present invention that receives and processes the leads as determined by the following sequencing stages. Deletions are "holes" (nucleotide deletions) in reads relative to the reference. This type of mismatch is referred to herein as a "d-type" mismatch. As with the classes "N" and "M", it is possible and appropriate to define a limit on the inaccuracy of the matching. The definition of the set of constraints for "Class I" is based on the same principles used for "Class M" and is shown in the last row of Table 1. In addition to the thresholds for each type of mismatch allowed for class I data, a further constraint is the number of mismatches "n", "s", "d", "i" and "c". , W (n, s, d, i, c). Such additional constraints make it possible to generate classes with an upper bound on matching inaccuracies according to any meaningful user-defined selective access criteria. For example, but not limited to, w (n, s, d, i, c) is allowed for leads belonging to (n + s + d + i + c) 1/5 or (n + s + d + i + c), or "class I". It can be any linear or non-linear expression that bounds the maximum inaccuracy level of the matching. Such a boundary is such that, beyond a simple threshold applied to each type of acceptable mismatch, this boundary is defined for any possible combination of the number of allowable mismatches in a "Class I" lead. It constitutes a very useful criterion for applying the desired selective data search to the compressed data when analyzing sequence reads for various purposes, to allow further boundaries to be set.
5. The fifth category is a mapping that is considered valid for each data class when referencing the reference sequence (ie, a set of matching rules that defines the upper limit of the maximum precision of the matching specified in Table 1). (Not met) include all leads that do not find. Such a sequence is referred to as “Unmapped” when referring to the reference sequence and is classified as belonging to “Class U”.
[Classification of lead pairs by matching rules]

The classification specified in the previous section is for a single sequence read.
For sequencing techniques (Illumina Inc.) that generate reads in pairs where the two reads are known to be separated by an unknown sequence of variable length, the entire pair is classified into a single data class. It is appropriate to consider doing so. A lead that is combined with another lead is called its "mate."

If both leads of a pair belong to the same class, the assignment of the entire pair to the class is obvious: the entire pair is assigned to the same class of any class (ie, P, N, M, I, U ). If two leads belong to different classes, but neither belong to "class U", the entire pair is assigned to the class with the highest priority defined according to the following equation:
P <N <M <I
Here, “class P” has the lowest priority and “class I” has the highest priority.

  If only one of the leads belongs to "class U" and its mate belongs to any of classes P, N, or M, the sixth class is defined as "class HM", which represents "half mapping".

  The definition of such a particular class of reads is based on the fact that it is used to try to determine gaps or unknown regions (also known as little-known unknown regions) present in the reference genome. . Such regions are reconstructed by mapping pairs at edges using pair reads that can map to known regions. The unmapped mate is used to create a so-called "contig" of the unknown area, as shown in FIG. Therefore, providing selective access to only these types of read pairs greatly reduces the associated computational burden and the use of modern solutions requires a large amount of data due to large data sets that need to be fully examined. Very efficient processing becomes possible.

  The following table shows the matching rules that are applied to leads to define the class of data to which each lead belongs. This rule is defined in the first five columns of the table regarding the presence or absence of mismatch types (n, s, d, i, c type mismatch). The sixth column provides the maximum threshold for each mismatch type and the rules for any functions f (n, s) and w (n, s, d, i, c) of the possible mismatch types. I do.

  Table 1. A set of mismatch types and restrictions that must be met in order for each sequence read to fall into the class of data defined in the present disclosure.

［マッチング精度の異なるサブクラスへのクラスＮ、Ｍ及びＩのシーケンスリードのマッチング規則のパーティション］ Table 1. A set of mismatch types and constraints that must be satisfied for each sequence read to be classified into the data classes defined in the present disclosure.

[Partitions of matching rules for class N, M and I sequence reads into subclasses with different matching accuracy]

The data classes of types N, M and I defined in the previous section can be further broken down into any number of distinct subclasses with different degrees of matching accuracy. Such an option is an important technical advantage in providing finer granularity and, as a result, provides more efficient selective access to each data class. By way of example and not limitation, the number of subclasses of class N k (subclass _{N 1,} · · ·, subclasses _{N k)} to decompose to the corresponding component _{MAXN 1, MAXN 2, ···,} MAXN (k-1) , MAXN _(k) must be defined, and under the conditions MAXN ₁ <MAXN ₂ <... <MAXN _(k−1) <MAXN, each lead is evaluated when each element of the vector is evaluated. Assign to the lowest ranked subclass that meets the restrictions specified in Table 1. This is illustrated in FIG. 29, where the data classification unit 291 includes classes P, N, M, I, U, HM encoders, and encoders for annotations and metadata. A class N encoder consists of a vector of thresholds from MAXN ₁ to MAXN _k 292 that generate k subclasses of N data (296).

For classes of type M and type I, the same principle applies by defining vectors with the same properties in MAXM and MAXTOT, respectively, with functions f (n, s) and w (n, s, d, I, Each vector component is used as a threshold to check if c) satisfies the constraint. As in the case of type N subclasses, the assignment is given to the lowest subclass for which the restriction is satisfied. The number of subclasses for each class type is independent, and any combination of subclasses is allowed. This is illustrated in Figure 29, class M encoder 293 and Class I encoder 294 are each configured MAXM _j from the threshold MAXM _1, and from MAXTOT ₁ vector of MAXTOT _h. The two encoders generate j subclasses of M data (297) and h subclasses of I data (298), respectively.

  If two leads of a pair fall into the same subclass, the pair belongs to the same subclass.

If the two reads of a pair are classified into different classes of subclasses, the pair belongs to a higher priority class subclass according to the following formula:
N <M <I
Here, N has the lowest priority and I has the highest priority.

If the two leads belong to different subclasses of any of the classes N, M, or I, the pair belongs to the subclass with the highest priority according to the following formula:
N ₁ <N ₂ <... <N _k
M ₁ <M ₂ <... M _j
I ₁ <I ₂ <... <I _h
Here, the highest index has the highest priority.
[Conversion of "external" reference sequence]

  Mismatches found in leads classified into classes N, M, I can be used to create a "deformation" reference that is used to more efficiently compress the representation of the lead.

A lead classified as belonging to class N, M or I (with respect to an "existing" (ie "external") reference sequence, denoted as RS ₀ ) will, according to the occurrence of an actual mismatch with the "translated" reference, it can be encoded with respect to "conversion" Reference sequence RS _1. For example, for read ^M _in belonging to class M (shown as the ith read of class M) that includes a mismatch with respect to reference sequence RS _n , read ^M _in = read ^P _{i (n + 1)} after “conversion” is A (Ref _n ) = Ref _{n + 1} . Here, A is converted from the reference sequence RS _n to the reference sequence _{RS n + 1.}

FIG. 19 shows that a read including a mismatch (belonging to class M) to reference sequence 1 (RS ₁ ) is completely corrected for reference sequence 2 (RS ₂ ) obtained from RS ₁ by correcting the base corresponding to the mismatch position. An example of a method of converting leads into matching leads is shown below. These remain classified and are coded along with other reads in the same data class access unit, but coding is done using only the descriptors and descriptor values required for class P reads. This transformation can be expressed as:

If the representation of transform A, which produces RS ₂ when applied to RS ₁ , plus the representation of lead pair RS ₂ corresponds to a lower entropy than the representation of the class M vs. TS ₁ lead,
Since a higher compression of data representation is achieved, it is advantageous to transmit a corresponding representation of the expression and lead pair RS ₂ conversion A.

The encoding of transform A for transmission in a compressed bitstream requires the definition of two additional descriptors, as defined in the table below.

  FIG. 26 shows an example of how reference transform is applied to reduce the number of mismatches coded in the mapped leads.

Note that in some cases, a transformation is applied to the reference.
A mismatch may occur in the representation of the lead that did not exist when referencing the reference before applying the transformation.
The type of mismatch can be changed, the lead contains A instead of G and all other leads contain C instead of G, but the mismatch remains in the same position.
-Different data classes and a subset of the data of each data class may refer to the same "transformed" reference sequence or a reference sequence obtained by applying a different transformation to the same existing reference sequence.

  FIG. 27 further illustrates that after the reference transformation has been applied and the lead has been represented using the “translated” reference, the lead is read from the appropriate descriptor set (eg, from class M using class P descriptors). Here is an example of how to change the type of coding from one data class to another class by encoding the read of This means that if, for example, the conversion changes all bases corresponding to the read mismatch of the bases actually present in the read, then the reads belonging to class M (the original non- "converted" reference sequence Occurs when a reference is virtually converted to a virtual lead of a class P virtual lead (when referring to a "converted" reference). Definitions of the set of descriptors used for each class of data are provided in the following sections.

FIG. 30 shows how different classes of data can recode the read using the same “transformed” reference R ₁ = A ₀ (R ₀ ) (300), or a different transformation A _N (301), A _M (302) and A _I (303) can be applied separately to each class of data.
[Definition of information required to express sequence read to descriptor block]

  Once the lead classification has been completed with the definition of the class, further processing, when represented as being mapped to a particular reference sequence, of a separate descriptor representing the remaining information that allows the reconstruction of the lead sequence Is to define a set. The data structure of these descriptors requires the storage of global parameters and metadata used by the decoding engine. These data are composed of a Genomic Dataset Header shown in the following table. A data set is defined as the set of elements of the code necessary to reconstruct genomic information relevant to a single genomic sequencing run and all of the following analyses. If the same genomic sample is sequenced twice in two separate runs, the resulting data is encoded as two separate data sets.

Table 1. Genome dataset header structure

A sequence read (ie, a DNA segment) that references a given reference sequence can be adequately represented by the following equation:
・ Start position on reference sequence (pos)
A flag to signal if the read has to be considered as a reverse complement to the reference (rcomp).
The distance (pair) to the mate pair in the case of pair read.
A read length (len) in the case of a sequencing technology that generates a variable read length. If the lead length is constant, the lead length associated with each lead can obviously be omitted and stored in the header of the main file.
・ For each mismatch:
-Position mismatch (nmis for class N, snpp for class M, indp for class I)
-Type of mismatch (not present in class N, snpt in class M, indt in class I)
Flags that indicate special characteristics of sequence reads, such as: Templates with multiple segments in sequencing.Each segment must be correctly aligned according to the aligner.Unmapped segments. Next segment • Signaling of first or last segment • Poor quality control • PCR or optical replication • Secondary alignment • Auxiliary alignment • When soft clipped nucleotide sequence is present (Class I indc)
A flag indicating the reference used for alignment and compression (eg, an “internal” reference of class U), if applicable (descriptor rtype).
-For class U, the descriptor indc uses an "internal" reference to identify the unmatched parts of the lead (usually edges) using a specified set of matching precision limits.
The ureads descriptor is used to directly encode reads that cannot be mapped to an available reference because they are pre-existing (ie, “external” reference genomes) or “internal” reference sequences.

  This classification produces a group of descriptors (descriptors) that can be used to unambiguously represent genomic sequence reads. The following table summarizes the descriptors required for each class of read aligned with an "external" (i.e., "existing") or "internal" (i.e., "built") reference. .

Table 2. A block of descriptors defined for each class of data

  Reads belonging to class P are characterized and completely re-created solely by position, reverse complement information, and the offset between the mate obtained by the sequencing technique that generates the mate pair, some flags, and the read length. Be composed.

  The next section describes in detail how these descriptors are defined for classes P, N, M, I, and class U is described in the following section.

Class HM applies only to lead pairs, a special case where one lead belongs to class P, N, M or I and the other lead belongs to class U.
[Position descriptor]

  In the position (pos) block, only the coded first read mapping position is stored as an absolute value on the reference sequence. All other location descriptors assume a value representing the difference to the previous location. Such modeling of the information source defined by the sequence of lead position descriptors is generally characterized by reduced entropy, especially for sequencing processes that produce high coverage results.

For example, FIG. 1 shows a method of describing the start position of the first alignment as the position “10000” on the reference sequence, and then describing the position of the second lead starting at the position 10180 as “1080”. At high coverage (> 50 ×), most of the position vector descriptors exhibit a high value of low values, such as 0 or 1, or other small integers. FIG. 1 shows how the positions of three read pairs are described in a pos block.
[Reverse complementary descriptor]

  Each read of a read pair generated by the sequencing technique can originate from any genome strands of the sequenced organic sample. However, only one of the two strands is used as a reference sequence. FIG. 2 shows a state in which one lead (lead 1) starts from one strand and the other lead (lead 2) starts from the other strand in the lead pair.

  When strand 1 is used as a reference sequence, lead 2 is encoded as the reverse complement strand of the corresponding fragment on strand 1. This is shown in FIG.

  In the case of coupled reads, there are four combinations of direct and reverse complement pairs. This is shown in FIG. The rcomp block encodes four possible combinations.

The same coding is used for the reverse complement information of the leads belonging to classes N, M, P, I. In order to allow selective access to different data classes, the reverse complement information of the leads belonging to the four classes is coded in different blocks as shown in Table 2.
[Pairing information descriptor]

The pair descriptor is stored in the pair block. Such a block stores a descriptor that encodes the information needed to reconstruct the original read pair when the applied sequencing technique generates a pair-wise read. At the time of this disclosure, most of the sequencing data is generated using techniques to generate paired reads, but not in all techniques. This is why the presence of this block is not necessary to reconstruct all the sequencing data information, if the genomic data sequencing technique considered does not generate the paired read information.
[Definition]

Mate pair: a lead associated with another lead in the lead pair (eg, in the example above, lead 2 is the mate pair of lead 1).
Pairing distance: from a position of the first read (pairing anchor, eg, the last nucleotide of the first read) to a second read (eg, the first nucleotide of the second read). The number of nucleotides placed on a reference sequence that separates a position.
-Most likely pairing distance (MPPD): This is the most likely pairing distance expressed in number of nucleotides.
Pairing distance location (PDD): PDD is a method of expressing the pairing distance by the number of leads that separate the leads from each mate present in a particular location descriptor block.
Most likely pairing distance location (MPPD): The most likely number of reads that separates leads from mate pairs present in a particular location descriptor block.
Pairing Error Location (PPE): Defined as MPPD or the difference between MPPD and the actual location of the mate.
A pairing anchor: the last nucleotide position of the first read in the pair, used as a reference to calculate the distance of the mate pair in terms of the number of read positions or the number of nucleotide positions.

  FIG. 5 shows a method of calculating a pairing distance between read pairs.

  The pair descriptor block is a pairing error calculated as the number of reads that are skipped to reach the mate pair of the first read of the pair, for a defined decoded pairing distance. Is a vector.

  FIG. 6 shows an example of how the pairing error is calculated by both the absolute value and the derivative vector (characterized by lower entropy for high coverage).

The same descriptor is used for the pairing information of the leads belonging to the classes N, M, P and I. To allow selective access to different data classes, as shown in FIGS. 8 (Class N), FIGS. 10, 12 and 14 (Class M), and FIGS. The pairing information of the leads belonging to the class is coded in different blocks.
[Pairing information for reads mapped on different reference sequences]

Mapping a first read to one reference sequence (eg, a first chromosome) and mapping a second read to another reference sequence (eg, a fourth chromosome) in the process of mapping the sequence reads to a reference sequence. Is not uncommon. In this case, the above pairing information needs to be integrated with additional information related to the reference sequence used to map one of the leads. This is achieved by coding:
A predetermined value (indicating that it is mapped to two different sequences (if the lead 1 or lead 2 is mapped onto a sequence that is not currently coded, it will show different values) flag).
2. A unique reference identifier that refers to the coded reference identifier in the main header structure, as shown in Table 1.
3. The third element contains mapping information for the reference identified at point 2 and represented as an offset to the last coded position.

  FIG. 7 shows an example of this scenario.

In FIG. 7, the genomic encoder signals this information by creating an additional descriptor in the paired block, since read 4 is not mapped onto the currently coded reference sequence. In the following example, the lead 4 of the pair 2 has the reference No. 4, but the currently coded reference is no. It is one. This information is coded using three components:
1) One special predefined value is coded as the pairing distance (in this case, 0xffffff).
2) The second descriptor provides the reference ID described in the main header (4 in this case).
3) The third element contains the mapping information of the associated reference (170).
[Class N read mismatch descriptor]

Class N includes all leads with only "n-type" mismatches,
At the location of the A, C, G or T base, the called base is found as N. All other bases in the read match perfectly with the reference sequence.

FIG. 8 shows the method:
The position of "N" in lead 1 is
The absolute position of lead 1 or
The derivative position relative to the previous "N" in the same read,
Coded as
The position of “N” on lead 2 is
-The absolute position of the length of lead 1 + lead 2, or
-Coded as the derivative position relative to the previous "N". In the nmis block, encoding of each read pair ends with a special "separator" symbol.
[Descriptors encoding substitutions (mismatches or SNPs), insertions, deletions]

  A substitution is defined as the presence of a different nucleotide base in a mapped read relative to that occurring at the same position in the reference sequence.

FIG. 9 shows an example of replacement in a mapped read pair. Each substitution is coded as a "position" (snpp block) and a "type" (snpt block). Depending on the statistical occurrence of substitutions, insertions or deletions, different source models of the relevant descriptors can be defined, and the symbols generated in the relevant blocks can be coded.
[Source Model 1: Position and Substitution as Type]
[Replacement position identifier]

The replacement position is calculated in the same manner as the value of the mmis block. That is,
In lead 1, the replacement is
As an absolute position of lead 1 or as a derivative position relative to a previous replacement of the same lead
Coded.
In lead 2, the replacement is
As the absolute position of the length of lead 2 + lead 1, or
Coded.

FIG. 10 shows how the permutation (if the symbol of the lead differs from the symbol of the reference sequence at the specified mapping position) is coded.
1. The position of the mismatch-with respect to the starting position of the read, or-with respect to the previous mismatch (derivative coding)
2. Mismatch types represented as codes calculated as shown in FIG.

In the snpp block, the coding of each read pair ends with a special "separator" symbol.
[Replacement descriptor]

  In the case of class M (and I, as described in the next section), a mismatch is determined from the actual symbol present in the reference by the corresponding substitution present in lead {A, C, G, T, N, Z}. The symbol is coded by an index (moving from right to left). For example, if the aligned read presents a C instead of a T located at the same position in the reference, the mismatch indicator is indicated as "4". The decoding process reads the encoded descriptor and moves the nucleotide at the specified position on the reference from left to right to obtain the decoded symbol. For example, “2” received for a position where G exists in the reference is decoded as “N”. FIG. 11 shows all possible substitutions and the respective coding symbols. Clearly different context-adaptive probability models can be assigned to each permutation index according to the statistical properties of each permutation type of each data class to minimize descriptor entropy.

  When employing the IUPAC ambiguity code, the permutation mechanism yields exactly the same result, but the permutation vector is S = ｛A, C, G, T, N, Z, M, R, W, S, Y, K, V, H, D, B}.

  FIG. 12 shows an example of coding permutations within a snpt block.

Some examples of permutation coding when an IUPAC ambiguity code is employed are shown in FIG. 13, and another example of a permutation index is shown in FIG.
[Encoding Insertion and Deletion]

  For class I, mismatches and deletions are coded by index from the actual symbol present in the reference to the corresponding replacement symbol present in the lead: {A, C, G, T, N, Z} ( Move from right to left). For example, if the aligned read shows a C instead of a T at the same position in the reference, the mismatch indicator will be "4". If the lead presents a deletion due to the presence of "A" in the reference, the coded symbol will be "5". The decoding process reads the coded descriptor, the nucleotide at a given position on the reference, and moves from left to right to search for the decoded symbol. For example, a “3” received for a location where a G is in the reference is decoded as a “Z”.

  Insertions are coded as 6, 7, 8, 9, 10 for the inserted A, C, G, T, N, respectively.

FIG. 15 shows an example of a method for encoding substitutions, insertions and deletions in class I read pairs. To support the entire set of IUPAC ambiguity codes, the permutation vector S = {A, C, G, T, N, Z} is defined as S = {A, C, G, T, N, Z, M, R, W, S, Y, K, V, H, D, B}. In this case, if the permutation vector has 16 elements, the insertion code must have different values, ie, 16, 17, 18, 19, 20. This mechanism is shown in FIG.
[Source Model 2: Permutation types and indels per block]

  For some data statistics, coding models different from those described in the previous section allow for the development of permutations and indels that generate low entropy sources. Such a coded model is an alternative to the techniques described above for mismatches only and mismatches and indels.

  In this case, one data block is defined for each replaceable symbol (IUPAC code excluding 5, IUPAC code with 16), and one block is defined for deletion and four blocks for insertion. For simplicity, the description focuses on the case where IUPAC code is not supported.

FIG. 17 shows how each block contains a single-type insertion or mismatch location. If there is no mismatch or insertion of that type in the coded read pair, a 0 is coded in the corresponding block. The header of each access unit includes a flag indicating the first block to be decoded so that the decoder can start decoding the blocks described in this section. In the example of FIG. 18, the first element to be decoded is position 2 of the C block. If a particular type of mismatch or indel does not exist in the read pair, a 0 is added to the corresponding block. On the decoding side, if the decoding pointer of each block points to a value of 0, the decoding process moves on to the next read pair.
[Encoding of additional signaling flags]

Each of the data classes (P, M, N, I) introduced above may require encoding of additional information regarding the nature of the encoded lead. This information may be, for example, related to a sequencing experiment (eg, indicating the likelihood of one read overlapping) or may represent some characteristic of the read mapping (eg, the first of a pair). Or second). In the context of the present invention, this information is coded in a separate block for each data class. The main advantage of such an approach is that this information can be selectively accessed only when needed and only in the required reference sequence area. Another example of using such a flag is as follows.
• Lead pairs • Leads mapped as appropriate pairs • Unmapped leads or mates • Leads from reverse strand or first / second of mate pairs • Not primary alignment • Platforms / vendors that failed to read PCR or optical duplication • Auxiliary alignment [construction of class U descriptors and “internal” references of unmapped “class U” and “class HM” reads]

  For leads belonging to class U or unmapped "HM class" paired leads, they satisfy a specified set of matching accuracy constraints to belong to any of classes P, N, M, or I " Since it cannot be mapped to an "external" reference sequence, one or more "internal" reference sequences are "built" and used for the compressed representation of reads belonging to these data classes.

For example, there are several ways to build a proper "internal" reference, rather than the following restrictions:
Splitting unmapped reads into clusters containing reads that share a common contiguous genomic sequence of at least the smallest size (signature). Each cluster can be uniquely identified by its signature, as shown in FIG.
Sorting leads in a meaningful order (eg, lexicographical order) and using the last N leads as an “internal” reference to N + 1 coding. This method is shown in FIG.
-On a subset of the leads of class U, so-called "so-called", so that all or related subsets of the leads belonging to that class can be aligned and coded according to the specified matching accuracy constraints or the new constraint set. Doing "de-novo assembly".

If the encoded lead can be mapped to an "internal" reference that meets the specified set of matching precision constraints, the information needed to reconstruct the lead after compression uses the following types of descriptors: Coded.
1. Start position of the matching portion of the internal reference from the viewpoint of the read number of the internal reference (pos block). This location can be coded as an absolute or derivative value relative to a previously coded lead.
2. Offset of the start position from the start point of the corresponding read of the internal reference (pair block). For example, if the lead length is constant, the actual position is pos * length + pair.
3. There may be mismatches coded as the location of the mismatch (snpp block) and the type (snpt block).
4. Portions of the lead that do not match (or match, but have a number of matches above a defined threshold) the internal reference (generally edges identified in pairs) are coded in indc blocks. As shown in FIG. 24, a padding operation can be performed on some edges of the internal reference used to reduce the entropy of mismatch coded in the indc block. The encoder can select an optimal padding strategy according to the statistical characteristics of the genomic data being processed. The padding strategies that can be selected are as follows.
a. No padding b. A fixed padding pattern selected according to the frequency of the currently coded data c. Variable padding patterns defined according to the statistical properties of the current context, defined for the last N coded leads, The specific kind of padding strategy may be signaled by a special value in the indc block header.
5. Flag (rtype block) indicating whether the read is coded internally self-generated, external or without reference.
6. Verbatim coded leads.

  FIG. 24 shows an example of such an encoding procedure.

  FIG. 25 shows an alternative encoding of unmapped leads on the internal reference, where the pos + pair descriptor has been replaced with a coded pos. In this case, pos represents the distance of the leftmost nucleotide position of read n to the leftmost nucleotide position of read n-1, relative to the position on the reference sequence.

  If the leads of class U are of variable length, an additional descriptor rlen is used to store the length of each lead.

  This coding approach can be extended to support N starting positions per lead so that the lead can be split into two or more reference positions. It encodes reads generated by sequencing technologies (such as Pacific Bioscience) that generate very long reads (50K + bases) that typically display the repetitive pattern generated by the sequencing methodology loop. Especially useful to do. The same approach can also be used to encode chimeric sequence reads, defined as reads with little or no overlap between two different parts of the genome.

The above approach is clearly applicable beyond the simple class U and can be applied to any block that contains a descriptor associated with a lead position (pos block).
[Alignment score descriptor]

  The mscore descriptor provides a score for each alignment. In the context of the present invention, used to represent per-read mapping / alignment scores generated by the Genome Sequence Read Aligner.

  The score is represented using an exponent and a mantissa. The number of bits used to represent the exponent and the mantissa are transferred as configuration parameters. By way of example, but not by way of limitation, Table 2 shows how this is defined in IEEE RFC 754 for an 11-bit exponent and a 52-bit mantissa.

The score for each alignment can be expressed as:
-1-bit code (S)
-11-bit exponent (E)
-53-bit mantissa (M)

Table 2. The alignment score can be expressed as a 64-bit double-precision floating-point value

Since the base (radix) used in the calculation of the score is 10, the result is as follows.
Score = −1 ^s × 10 ^E × M
[Group of Leads]

Various types of sequence reads can be generated during the sequencing process. By way of example, and not limitation, types may be associated with different sequenced samples, different experiments, and different configurations of sequencing equipment. According to the disclosed invention, this information is stored after sequencing and alignment by a dedicated descriptor named rgroup. rgroup is a label associated with each coded lead, and allows the decoding device to divide the decoded leads into groups after decoding.
[Multiple alignment descriptor]

The following descriptors are specified to support multiple alignment: If there is a spliced lead, the present invention defines spline_reads_flag as a global flag that is set to one.
[Mmap descriptor]

The mmap descriptor is used to notify how many positions of the left end of a lead or a pair of leads are aligned. Genome records containing multiple alignments are associated with one multi-byte mmap descriptor. The first two bytes of the mmap descriptor may be a single segment (if there is no splice in the encoded data set) or, alternatively, all segments where the read was spliced for some possible alignments An unsigned integer N that refers to the read as (if a splice exists in the data set). The value of N indicates how many pos descriptor values are coded for the template for this record. As described below, one or more unsigned integer M _i after N is followed.
[Strand property of multiple alignment]

The rcomp descriptor described in the present invention is used to specify the strandedness of each read alignment using the syntax specified in the present invention.
[Multiple alignment score]

In the case of multiple alignment, one mscore specified in the present invention is assigned to each alignment.
[Multiple alignment without splice]

  If there is no splice in the access unit, the setting of spline_reads_flag is released.

In paired-end sequencing, the mmap descriptor is a 16-bit unsigned integer N followed by i, assuming i takes a value from 1 to the number of the completely first (here, left-most) read alignment alignment. It consists of an integer M _i unsigned one or more 8-bit. For the first read alignment, whether spliced or not, Mi is the number of segments used to align the second read (in this case, if there is no splice, this is the number of alignments). Equals) and is used to signal how many pair descriptor values are coded for the alignment of the first read.

が使用される。 The value of M _i, to represent the number of the second lead of the alignment, the following equation

Is used.

The special value of M _i (= 0) is that the i-th alignment of the left-most lead is the k-th alignment of the left-most lead with k <i (when no new alignment matching the above equation is detected). Indicates that it is paired with the alignment of the right-most lead that has already been paired.

For example, in the simplest case:
1 If there is a single alignment for the left end lead and two alternative alignments for the right end, N is 1 and M1 is 2.
2 two alternative alignment is detected in the left end of the lead, if only detected one at the right end of the lead, N is the 2 next, M ₂ is zero.

When M _i is 0, the associated value pairs must be linked to the existing second lead alignment; otherwise the syntax error occurs, is regarded as the alignment is broken.

Example: As mentioned above, if the first lead has two mapping positions and the second lead has one mapping position, N is 2, M ₁ is 1, and M ₂ is 0. Subsequent to this, if another alternative secondary mapping to the entire template is performed, N is 3 and M ₃ is 1.

39, N in the case of no multiple alignment of splice and error, P, indicates the meaning of _{M i,} the reference source is not found, encode multiple alignment information using pos, pair, and a mmap descriptor Here's how.

のアライメントを有し、
・ 左端のリードの第ｉ番目のアライメントが、左端のリードの第ｋ番目（ｋ＜ｉ）のアライメントと既にペアになっている右端のリードのアライメントと、ペアになっている場合のＭ_ｉのいくつかの値は＝０になることがあり、
・ ペア記述子の１つの予め定められた値は、他のＡＵの範囲に属するアライメントの信号に存在することができる。それが存在する場合は、常に、現在のレコードに対する第１のｐａｉｒ記述子になる。
［スプライスを使用したマルチプルアライメント］ For 40, it is as follows:
・ The lead on the right end

Has the alignment of
- the left end of the i-th alignment of the lead, the k-th of the left end of the lead (k <i) the alignment of the already of the right end of the lead in the pair alignment, of if it is to a pair of M _i Some values can be = 0,
-One predetermined value of the pair descriptor can be present in the signal of the alignment belonging to the range of another AU. If it exists, it will always be the first pair descriptor for the current record.
[Multiple alignment using splice]

  If the dataset is coded with spliced reads, the msar descriptor can be used to represent splice length and strandedness.

  After decoding the mmap and msar descriptors, the decoder knows the number of leads or lead pairs coded to represent multiple mappings and the number of segments that make up the mapping of each lead or lead pair. This is shown in FIG. 41 and FIG.

のスプライスを有し、ここでＭ_ｉは、左端のリードの第ｉ番目のアライメント（１≦ｉ≦Ｎ_１）とペアで関連付けられた右端のリードのスプライスの数である。つまり、Ｐは右端のリードのスプライスの数を表し、ｍｍａｐ記述子の最初の値に続くＮ値を使用して計算される。
・ Ｎ_１及びＮ_２は、第１及び第２のリードのアラインメントの数を表し、ｍｓａｒ記述子のＮ＋Ｐ値を使用して計算される。 Referring to FIG. 41, the following applies:
The-left end of the lead, having _{N 1} alignment with the N splice _{(N 1} ≦ N).
N represents the number of splices present in all alignments of the leftmost lead, coded as the first value of the mmap descriptor.
・ The right end lead

Where M _i is the number of splices in the rightmost lead paired with the ith alignment of the leftmost lead (1 ≦ i ≦ N ₁ ). That is, P represents the number of splices of the rightmost lead and is calculated using the N value following the first value of the mmap descriptor.
· N ₁ and _{N 2} represents the number of alignment of the first and second leads, are calculated using the N + P value of msar descriptor.

スプライス、ｔ_ｊ １≦ｊ≦Ｐ、及びＮ_２（Ｎ_２≦Ｐ）アライメントを有する。
・ ｐａｉｒ記述子の数は、ＮＰ＝Ｍａｘ（Ｎ１，Ｐ）＋Ｍ_０として計算され、ここで
・ Ｍ０は値が０のＭｉの数であり、
・ ＮＰは、1つの特別なｐａｉｒ記述子が他のＡＵにアラインメントが存在することを示す場合に１だけ増加する必要がある。
［アライメントスコア］ Referring to FIG. 42, the following applies:
- the left end has a _{N 1} alignment with the N splice _{(N 1} ≦ N). If N ₁ = N AND N ₂ = P, there is no splice.
・ The right end lead

With splice, t _j 1 ≦ j ≦ P, and N ₂ (N ₂ ≦ P) alignment.
· Pair number of descriptors are calculated as NP = Max (N1, P) + M 0, where · M0 is the number of values is 0 Mi,
NP needs to be incremented by 1 if one special pair descriptor indicates that there is an alignment in another AU.
[Alignment score]

The mscore descriptor allows notification of the alignment mapping score. In single-ended sequencing, having N ₁ value for each template; in paired end sequencing has values for each alignment of the entire template (of the first lead different alignment numbers + the second lead further Number of alignments, ie, if M _i −1> 0)
Number of scores = MAX (N ₁ , N ₂ ) + M ₀
Here, M0 indicates the total number of M _i = 0.

In the present invention, a plurality of score values can be associated with each alignment. The number of alignments is signaled by the configuration parameter as_depth.
[Descriptor for multiple alignment without splice]

［スプライスを使用したマルチプルアライメントの記述子］ Table 3. Determining the number of descriptors needed to represent multiple alignments within a single genome record for spliceless multiple alignments

[Descriptor for multiple alignment using splice]

  Table 4 shows the determination of the number of descriptors needed to represent multiple alignments in one genomic record, for multiple alignments with splices.

［異なるシーケンス上のマルチプルアラインメント］ Table 4. Determining the number of descriptors needed to represent multiple alignments in one genomic record in multiple alignments with splices

[Multiple alignment on different sequences]

  The alignment process may find an alternative mapping to a reference sequence different from the reference sequence where the primary mapping is located.

  For uniquely aligned read pairs, for example, when there is a chimeric sequence with a mate on another chromosome, the pair descriptor must be used to represent the absolute read position. The pair descriptor must be used to signal the reference and the location of the next record that contains further alignment to the same template. The last record (eg, the third if the alternate mapping is coded with three different AUs) includes the reference and the location of the first record.

If one or more alignments of the leftmost lead of the pair are on a different reference sequence than the reference sequence associated with the currently coded AU, a predetermined value is used for the pair descriptor. The predetermined value is followed by the reference sequence identifier and the position of the leftmost alignment in everything contained in the next AU (ie, the first decoded value of the pos descriptor of the record).
[Multiple alignment including uninserted, deleted and unmapped parts]

If the alternative secondary mapping does not preserve the continuity of the reference region to which the sequence is aligned, the actual sequence (and descriptors associated with mismatches such as substitutions or indels) will be It may not be possible to reconstruct the exact mapping generated by the aligner, since it is only coded for this. The msar descriptor is used to indicate how the secondary alignment is mapped to the reference sequence when indels and / or soft clips are included. If msar is represented by the special symbol “*” of the secondary alignment, the decoder reconstructs the secondary alignment from the mapping positions of the primary alignment and the secondary alignment.
[Msar descriptor]

  The msar (Multiple Segments Alignment Record) descriptor supports spliced leads and alternative secondary alignments including indels or soft clips.

msar aims to communicate the following information:
• the length of the mapped segment; • the serialization of different mappings of the secondary alignment and / or spliced reads (ie, the presence of inserted, deleted or clipped bases).

  msar uses the extended CIGAR string syntax described below and the additional symbols described in Table 5.

［拡張シガー構文］ Table 5. Special symbols used in msar descriptors in addition to the syntax described in Table 6

[Extended cigar syntax]

  This section specifies the extended CIGAR (E-CIGAR) syntax for associating strings with information about sequences and associated mismatches, indels, clipped bases, multiple alignments, and spliced reads.

  Table 6 lists the editing operations described in the present invention.

［ソースモデル、エントロピーコード化及びコード化モード］ Table 6. MPEG-GE-CIGAR string syntax

[Source model, entropy coding and coding mode]

  For each data class, subclass and associated descriptor block of the genomic data structure disclosed in the present invention, different coding algorithms are employed according to the particular characteristics of the data or metadata obtained by each block and its statistical properties. May be. The "coding algorithm" needs to be intended as an association between a particular "source model" of the descriptor block and a particular "entropy coder". A particular "source model" can be specified and selected to obtain the most efficient coding of the data with respect to minimizing source entropy. The choice of entropy coder can be driven by coding efficiency considerations and / or probability distribution features and related implementation issues. Each selection of a particular "coding algorithm", also referred to as a "coding mode", can be applied to the entire "descriptor block" associated with the data class or subclass of the entire dataset, or divided into access units. A different "coding mode" can be applied to each part of the descriptor.

Each "source model" associated with a coding mode is characterized as follows:
The definition of descriptors originating from each source (ie, the set of descriptors used to represent the class of data, such as read position, read pairing information, mismatch to reference sequence, etc. as defined in Table 2) ).
• Definition of the associated probability model.
• Definition of the relevant entropy coding.
[Another advantage]

  Classification of sequence data into defined data classes and subclasses is characterized by modeling the sequence of descriptors by a single, discrete data source (eg, distance, location, etc.), lower information source entropy Enables efficient coding mode implementation using

  Another advantage of the present invention is that only a subset of the type of data of interest can be accessed. For example, one of the most important applications in genomics is to find differences in genomic samples relative to a reference (SNV) or population (SNP). Today, such analysis requires the processing of a complete sequence read, but by employing the data representation disclosed by the present invention, mismatches are already separated into only one to three data classes. (Depending on interest in also considering mismatches between "n-type" and "i-type").

  A further advantage is that when a new reference sequence is published, or when remapping is performed on data that has already been mapped to obtain a new alignment (eg, using a different mapping algorithm), certain It is possible to perform efficient transcoding from data and metadata compressed with reference to an “external” reference sequence to another different “external” reference sequence.

FIG. 20 shows a coding device 207 based on the principles of the present invention. The coding device 207 receives, for example, raw sequence data 209 generated by the genome sequencing device 200 as an input. Genome sequencing device 200 is well known in the art, such as the Illumina HiSeq 2500 (Illumina HiSeq 2500) or Thermo-Fisher Ion Torrent device. The raw sequence data 209 is supplied to the aligner unit 201 to prepare a sequence to be encoded by aligning the read with the reference sequence 2020. Alternatively, using a dedicated module 202, using different strategies as described in the sections "Building Internal References for Unmapped Leads of Class U" and "HM Classes" herein. , A reference sequence can be generated from the available leads. After being processed by the reference generator 202, the reads may be mapped onto the resulting longer sequence. The aligned sequence is then classified by the data classification module 204. Next, a further step of reference conversion is applied to the reference to reduce the entropy of the data generated by the data classification unit 204. This means that the external reference 2020 is processed by the reference conversion unit 2019 that generates the converted data class 2018 and the reference conversion descriptor 2021. Next, the converted data class 2018 is supplied to the block encoders 205 to 207 together with the reference conversion descriptor 2021.
The genomic block 2011 is then provided to arithmetic encoders 2012-2014 that encode the block according to the statistical properties of the data or metadata carried by the block. The result is a genome stream 2015.

  FIG. 21 shows a decoding device 218 according to the principles of the present disclosure. The decoding device 218 receives the multiplexed genomic bitstream 2011 from a network or a storage element. The multiplexed genomic bit stream 2110 is provided to a demultiplexer 210 to generate individual streams 211, which are then provided to entropy decoders 212-214, where the genomic blocks 215 and a reference transform descriptor 2112 is generated. The extracted genomic blocks are supplied to block decoders 216 to 217, the blocks are further decoded into data classes, and a reference conversion descriptor is supplied to a reference conversion unit 2113. The class decoder 219 further processes the genomic descriptor 2111 and the transformed reference 2114 and merges the results to generate an uncompressed read of the sequence, which is further converted to a format known in the art, such as a text file. Alternatively, it can be stored in a zip compressed file, or a FASTQ or SAM / BAM file.

  The class decoder 219 converts the original genomic sequence by utilizing information about the original reference sequence carried by one or more genomic streams, and the reference transform descriptor 2112 carried in the encoded bitstream. Can be rebuilt. If reference sequences are not transferred by the genomic stream, they must be available on the decoding side and accessible by the class decoder.

  The techniques of the present invention disclosed herein can be implemented in hardware, software, firmware, or any combination thereof. If implemented in software, they may be stored on a computer medium and executed by a hardware processing unit. A hardware processing unit may include one or more processors, digital signal processors, general purpose microprocessors, application specific integrated circuits or other discrete logic circuits.

The techniques of this disclosure may be implemented with various devices or devices, including mobile phones, desktop computers, servers, tablets, and similar devices.
[File format: Selective access to genomic data area using master index table]

To support selective access to specific regions of aligned data, the data structures described herein implement an indexing tool called a master index table (MIT). This is a multi-dimensional array containing the locations where a particular read is mapped to the associated reference sequence.
The value included in the MIT is the mapping position of the first read in each pos block so that non-sequential access to each access unit is supported. The MIT includes a section for each class of data (P, N, M, I, U, and HM) and a reference sequence. The MIT is included in the genomic data set header of the encoded data. FIG. 21 shows the structure of a Genomic Dataset Header, FIG. 32 shows a general visual representation of MIT, and FIG. 33 shows an example of MIT for class P of coded reads.

  The values contained in the MIT shown in FIG. 33 are used to directly access the region of interest (and the corresponding AU) in the compression domain.

  For example, with reference to FIG. 33, if it is necessary to access an area included between locations 150,000 and 250,000 on reference 2, the decoding application skips to the MIT second reference and k1 Search for two values k1 and k2 such that <150,000 and k2> 250,000. Here, k1 and k2 are two indexes read from the MIT. In the example of FIG. 33, this is the second and third position of the second vector of the MIT. These return values are used by the decoding application to get the appropriate data location from the pos block's local index table, as described in the next section.

Along with pointers to blocks containing data belonging to the four classes of genomic data described above, the MIT can be used as an index for additional metadata and / or annotations added to the genomic data during its life cycle. it can.
[Local index table]

  At the beginning of each genome data block, a data structure called a local header is attached. The local header includes a unique identifier of the block, a vector of access unit counters for each reference sequence, a local index table (LIT), and optionally, block-specific metadata. LIT is a vector of pointers to the physical locations of the data belonging to each access unit in the block payload. FIG. 34 illustrates a typical block header and payload where the LIT is used to access a particular area of coded data in a non-sequential manner.

In the previous example, the reference sequence No. To access regions 150,000 to 250,000 of reads aligned to 2, the decoding application searched locations 3 and 4 from the MIT. These values are used by the decoding process to access the third and fourth elements of the corresponding section of the LIT. In the example shown in FIG. 35, the LIT index related to the AU related to reference 1 (5 in the example) is skipped using the total access unit (Total Access Units) counter included in the block header. Therefore, an index containing the physical location of the requested AU in the encoded stream is calculated as follows:
Position of data block belonging to requested AU = data block belonging to AU of reference 1 to be skipped + position searched using MIT First block position: 5 + 3 = 8
Last block position: 5 + 4 = 9

  A block of data obtained using an indexing mechanism called a local index table is part of the requested access unit.

  FIG. 26 shows how blocks included in the MIT table correspond to LIT blocks for each class or subclass of data.

  FIG. 37 shows how data blocks retrieved using MIT and LIT constitute one or more access units as defined in the next section.

In one embodiment of the present invention, the LIT can be integrated as a substructure of the MIT. The advantage of such an approach is the speed of accessing indexed data in the case of sequential parsing of the compressed file. If the LIT is integrated into the MIT of the file header, the decoding device need only analyze a small portion of the data to retrieve the requested compression information for selective access. Another advantage is that, when streaming over a network, the index information contained in the MIT and LIT is delivered in the first data block, thus allowing the receiving device to sort and select before the entire data transfer is completed. It is clear to a person skilled in the art to be able to perform operations such as dynamic access [access unit].

  Genomic data, categorized by data class and structured by compressed or uncompressed blocks, is organized into different access units.

  Genome access units (AUs) can be used to reconstruct nucleotide sequences and / or associated metadata (generally or uncompressed) and / or sequences of DNA / RNA (eg, a virtual reference). And / or a section of annotation data generated by a genome sequencing device and / or a genome processing device or analysis application. FIG. 37 shows an example of the access unit.

  An access unit is a block of data that can be decoded independently of other access units, using only globally available data (eg, a decoder configuration) or using information contained in other access units. .

Access units are distinguished as follows:
Characterize the type, the nature of the genomic data, the datasets they hold, and the way to access them,
Order, which provides a unique order for access units belonging to the same type.

All types of access units can be further divided into different “categories”.

The following is a non-exhaustive list of definitions of various types of genomic access units:
1) Type 0 access units need not refer to information from other access units that are accessed or decoded and accessed. All the information transmitted by the data or data sets they contain is read and processed independently by the decoding device or processing application.
2) Type 1 access units include data referencing data transmitted by type 0 access units. Reading or decoding and processing of data contained in Type 1 access units requires access to one or more Type 0 access units. Type 1 access units encode genomic data associated with “class P” sequence reads.
3) Type 2 access units include data referencing data transmitted by type 0 access units. Reading or decoding and processing of the data contained in the type 2 access units requires access to one or more type 0 access units. Type 2 access units encode genomic data associated with "class N" sequence reads.
4) Type 3 access units include data referencing data transmitted by type 0 access units. Reading or decoding and processing of data contained in Type 3 access units requires access to one or more Type 0 access units. Type 3 access units encode genomic data associated with “class M” sequence reads.
5) Type 4 access units include data referencing data transmitted by type 0 access units. Reading or decoding and processing of data contained in a type 4 access unit requires access to one or more type 0 access units. Type 4 access units encode genomic data associated with “class I” sequence reads.
6) Type 5 access units contain leads that cannot be mapped to available reference sequences ("Class U") or the like and are coded using internally built reference sequences. The type 5 access unit includes data that refers to data transmitted by the type 0 access unit. Reading or decoding and processing of data contained in a type 5 access unit requires access to one or more type 0 access units.
7) The type 6 access unit includes a read pair. One of the leads belongs to any of the classes P, N, M, and I, and the other read has an available reference sequence (“HM class”). ) Cannot be mapped. The type 6 access unit includes data referencing data transmitted by the type 0 access unit. Reading or decoding and processing of data contained in a type 6 access unit requires access to one or more type 0 access units.
8) Type 7 access units include metadata (eg, quality scores) and / or annotation data associated with the data or data set contained in the type 1 access units. Type 7 access units may be classified and labeled into different blocks.
9) Type 8 access units include data or data sets classified as annotation data. Type 8 access units may be classified and labeled on a block basis.
10) Additional access units can extend the structures and mechanisms described herein. As an example, but not by way of limitation, the results of genomic variant calling, structural and functional analysis can be encoded into a new type of access unit. The data organization in the access unit described herein does not preclude any kind of data being encapsulated in the access unit, which is a completely transparent mechanism with respect to the nature of the coded data.

  Type 0 access units are ordered (eg, numbered) and need not be stored and / or transmitted in an ordered manner (technical advantage: parallel processing / parallel streaming, multiplexing).

  Access units of types 1, 2, 3, 4, 5 and 6 do not need to be ordered and need not be stored and / or transmitted in an ordered manner (technical advantage: parallel processing / parallel streaming).

  FIG. 37 shows how an access unit is made up of one or more blocks of the same kind of data as a header. Each block can be composed of one or more blocks. Each block contains a number of packets, where the packets are a structured sequence of descriptors introduced above to represent, for example, lead positions, pairing information, reverse complement information, mismatch positions and types, etc. is there.

  Each access unit can have a different number of packets per block, but within an access unit, all blocks have the same number of packets.

Each data packet can be identified by a combination of three identifiers XYZ:
X indicates the access unit to which it belongs,
Y indicates the block to which it belongs (ie the type of data to be encapsulated),
-Z is an identifier indicating the packet order for other packets in the same block.

  FIG. 38 shows an example of an access unit and a packet label. Here, AU_T_N is a type T access unit having an identifier N, which may or may not imply the concept of the order depending on the type of the access unit. The identifier is used to uniquely associate one type of access unit with another type of access unit that is required to fully decrypt the transferred genomic data.

All types of access units are further classified and displayed in different "categories" according to different sequencing processes. For example, without limitation, classification and display can be performed in the following cases.
1. Sequencing the same organism at different times (access units contain genomic information with "temporary" meaning)
2. Sequencing organic samples of different properties of the same organism (samples of human skin, blood, hair, etc.) These are access units meaning "biological".

Claims (66)

A method of encoding genomic sequence data comprising a nucleotide sequence read, said method comprising:
Aligning said lead with one or more reference sequences, thereby creating an alignment lead;
Classifying said alignment reads using said one or more reference sequences according to a specified matching rule, thereby creating a class of alignment reads;
Encoding the classified alignment reads as multiple blocks of descriptors,
Encoding the classified alignment reads as multiple blocks of the descriptor includes selecting the descriptor according to the class of the alignment reads,
Structuring the block of descriptors with header information, thereby creating a continuous access unit. Further comprising classifying the leads that do not satisfy the specified matching rule into a class of unmapped leads,
Constructing a set of reference sequences using at least some of the unmapped reads;
Aligning the unmapped class of reads with the constructed set of reference sequences,
Encoding the classified alignment reads as multiple blocks of descriptors,
Encoding the set of constructed reference sequences;
Constructing the block of descriptors and the coded reference sequence with header information, thereby creating a continuous access unit;
The coding method according to claim 1. The classification includes classifying a genomic read with no mismatch in the reference sequence as a first `` Class P '' if there is no mismatch in the mapped read with respect to the reference sequence used for mapping.
The method according to claim 2. The classification is based on the genomic read if the sequencing device cannot call any "bases" and a mismatch is found only at positions where the number of mismatches in each read does not exceed a predetermined threshold. Further classifying as a second "class N".
The method of claim 3. The classification further comprises identifying the genomic read as a third "class M" if a mismatch is found at a position where the sequencing device was unable to call any "bases", wherein "n-type" ), And / or referred to as “bases” that are different from the reference sequence and are termed “s-type” mismatches, and the number of mismatches is greater than the “n-type” mismatch, the “s-type” mismatch. A predetermined threshold is not exceeded for the number of "type" mismatches, said threshold being given by a function (f (n, s)),
The method according to claim 4. The classification further includes identifying the genomic read as a fourth "Class I" if a mismatch of the same type as the "Class M" is likely to occur, wherein at least one mismatch type: " Inserts ("i-type"), "deletions"("d-type"), and soft clips ("c-type"), where the number of mismatches of each type is determined by a corresponding predetermined threshold And the threshold is given by a function (w (n, s, i, d, c)),
The method of claim 5. The classifying further includes identifying the genomic read as a fifth "class U" as including all reads that do not find any of the classes P, N, M, I.
The method of claim 6. Encoded genomic sequence reads are paired,
The method according to claim 7. The classification may be such that the classification includes genomic reads as a sixth read, with one read belonging to class P, N, M or I and the other read including all read pairs belonging to "class U". Further comprising identifying as "Class HM".
The method according to claim 8. Identifying whether the leads of the two mates are classified into the same class (P, N, M, I, U respectively) and assigning the pair to the same identified class;
Identify whether the two mate's leads are classified into different classes, and if neither belongs to "Class U", assign the paired leads to the highest priority class according to:
P <N <M <I
Here, "class P" has the lowest priority, "class I" has the highest priority,
Identifying whether only one of the two mate's leads is classified as belonging to "Class U" and classifying the paired leads as belonging to a "Class HM" sequence.
An encoding method according to claim 9. Each class of leads N, M, and I has the number of “n-type” mismatches (292), function f (n, s) (293), and function w (n, s, i, d, c) (294) Is further divided into two or more subclasses (296, 297, 298) according to the threshold vector (292, 293, 294) defined for each class N, M, I, respectively.
The method according to claim 11. Identifying whether the leads of the two mates are classified in the same subclass, assigning the pair to the same subclass,
Identifying whether the leads of the two mates are classified into different classes of subclasses and assigning the pair to the subclass belonging to the higher priority class according to the following formula:
N <M <I
Here, N has the lowest priority, I has the highest priority,
Identify whether the two mates' leads are classified into the same class, which class is N, M, or I, but have different subclasses, and assign the pair to the highest priority according to the following formula: Assign to higher subclass,
N ₁ <N ₂ <... <N _k
M ₁ <M ₂ <... M _j
I ₁ <I ₂ <... <I _h
Where the highest index has the highest priority,
The method according to claim 11. Information about the mapping location of each lead is coded by a "pos" descriptor block,
The method according to claim 12. Information about the strandiness of each read (ie, the sequence from which the DNA strand reads are derived) is encoded by the rcomp descriptor block,
The method according to claim 13. The pairing information of the paired-end read is encoded by a “pair” descriptor block.
The method according to claim 14. Additional alignment information, whether the reads are mapped in the correct pair, platform / vendor quality check failure, PCR or optical replication, or ancillary alignments, is included in the "flags" descriptor. Coded by blocks,
The method according to claim 15. Information about the unknown base is encoded by a "mmis" descriptor block,
The method of claim 16. Information about the location of the replacement is coded by an "snpp" descriptor block,
The method according to claim 17. Information about the type of substitution is coded by a particular "snpt" descriptor block,
The method according to claim 18. Information about the location of the mismatch, substitution, insertion or deletion is coded by an "indp" descriptor block,
The method according to claim 19. Information about the type of substitution, insertion, or deletion mismatch is encoded by an "indt" descriptor block;
The method according to claim 20. Information about the clipped bases of the mapped read is encoded by an "indc" descriptor block,
A method according to claim 21. Information about unmapped leads is coded by an “ureads” descriptor block,
23. The method according to claim 22. Information about the type of reference sequence used for encoding is encoded by an "rtype" descriptor block;
A method according to claim 23. Information about the multiple alignment of the mapped read is encoded by a "mmmap" descriptor block;
A method according to claim 24. Information about the spliced alignment and multiple alignment of the same read is encoded by an "msar" descriptor block and an "mmp" descriptor block;
A method according to claim 25. Information about the alignment score of the read is encoded by an “mscore” descriptor block;
The method of claim 26. Information about the group to which the lead belongs is coded by an “rgroup” descriptor block,
A method according to claim 27. A class P access unit is constructed using a block of descriptors of type "pos", "rcomp", and "flags".
29. The method according to claim 28. The access unit of class P encodes pair-end pairing information using a block of "pair"descriptors;
A method according to claim 29. Class N access units are constructed using the same block of class P access unit descriptors, in addition to using the "nmis" descriptor block for information on the location of unknown bases,
31. The method according to claim 30. A class M access unit is constructed using the same block of class P access unit descriptors in addition to the "snpp" and "snpt" descriptor blocks for information about the location and type of permutation.
31. The method according to claim 30. The class I access unit is the same as the class P access unit descriptor, in addition to the "indp", "indt", and "indc" descriptor blocks for information about substitutions, insertions, deletions, and the location and type of clip bases. Built using blocks,
31. The method according to claim 30. A class HM access unit is constructed using the same block of descriptors of the class I access unit for the mapped read and the block of the "ureads" descriptor for the unmapped read,
A method according to claim 33. Information about the multiple alignment is conveyed using blocks of “mmmap” and “msar” descriptors,
A method according to claim 33. Information about the spliced alignment
・ Symbol for displaying matching base =
・ Symbol for displaying insertion +
・ Symbol for displaying deletion ・ Symbol for displaying forward strand splice /
Symbol% for indicating the strand of the reverse strand
・ Symbol for indicating non-directional splice *
A text character from the IUPAC code for the DNA to indicate the substitution; a symbol (n) to display n soft-clip bases, where n is an integer; a symbol to display n hard-clipped bases. [N], where n is transmitted in an extended cigar string containing
A method according to claim 35. The block of descriptors includes a "master index table" including one section for each class and subclass of aligned reads, the sections being coded in both the "master index table" and the access unit. Including the mapping location on the one or more reference sequences of the first read of each access unit of each class or subclass of data.
37. The method of claim 36. The block of the descriptor further includes the type of reference used (existing or constructed) and information about the segment of the lead that is not mapped to the reference sequence;
38. The method of claim 37. The reference sequence is first transformed to a different reference sequence by applying substitution, insertion, deletion, and clipping, and the encoding of the classified alignment read as multiple blocks of descriptors is transformed. Referencing the reference sequence
39. The method according to claim 38. The same transformation is applied to the reference sequence used for all classes of data,
40. The method of claim 39. Different transformations are applied to the reference sequence used for each class of data;
41. The method according to claim 40. The transformation of the reference sequence is coded as a block of descriptors and structured with header information, thereby creating a continuous access unit;
42. The method of claim 41. The associated reference sequence transformation as a multiplex of the coding and descriptor blocks of the classified alignment read includes associating with a particular descriptor block and a particular source model.
43. The method according to claim 42. The entropy coder is one of a context adaptive arithmetic coder, a variable length coder, and a Golomb coder.
A method according to claim 43. A method of decoding encoded genomic data, said method comprising:
Analyzing the access unit containing the encoded genomic data to extract a multiplexed block of descriptors using the header information;
Decoding the multiplexed block of descriptors to extract leads according to specific matching rules that define a classification for one or more reference sequences. Further comprising decoding a master index table including one section for each class of associated associated mapping locations and leads;
The decoding method according to claim 45. The type of reference used: further including decoding information related to existing, transformed, or constructed;
47. The decoding method according to claim 46. Decoding the information associated with one or more transforms applied to the existing reference sequence.
A decoding method according to claim 47. The block of descriptors is entropy decoded,
49. The decoding method according to claim 48. Class P leads are obtained by decoding blocks of each type of descriptor, "pos", "rcomp", "flags", and "rlen",
Class N reads are obtained by decoding blocks of each type of descriptor, "pos", "rcomp", "flags", "rlen", "nmis",
Class M leads are obtained by decoding blocks of each type of descriptor, "pos", "rcomp", "flags", "rlen", "snpp", and "snpt",
Class I leads are obtained by decoding blocks of descriptors of each type: "pos", "rcomp", "flags", "rlen", "indp", "indt", and "indc";
The class U lead is a block of descriptors of each type of "pos", "rcomp", "flags", "rlen", "snpp", "snpt", "indc", "ureads", and "rtype". Obtained by decoding
50. The method according to claim 49. Classes P, N, M, and I are also obtained by decoding a block of "pair" descriptors,
The class HM is also obtained by decoding a block of descriptors "pos", "rcomp", "flags", "rlen", "indp", "indt", "indc", and "ureads". ,
A decoding method according to claim 50. Genome sequence data 209, a genomic encoder (210) for compressing said genomic sequence data 209 including nucleotide sequence reads, wherein the method comprises:
The genome encoder (210)
An aligner unit (201) configured to align the lead with one or more reference sequences, thereby creating an alignment lead;
A constructed reference generation unit (202) configured to generate the constructed reference sequence;
Data configured to classify said alignment reads according to a particular matching rule using one or more existing or constructed reference sequences, thereby creating a class of alignment reads (208). A classification unit (204);
One or more block coding units (205-207) configured to code the classified alignment read as a block of descriptors by selecting the descriptor according to the classified alignment read; ,
A multiplexer (2016) for multiplexing the compressed genomic data and metadata. A reference sequence conversion unit (2019) configured to convert the existing reference and data class (208) to the converted data class (2018);
A genome encoder according to claim 52. Said data classification unit (204) includes encoders of data classes N, M and I composed of a vector of thresholds producing subclasses of data classes N, M and I;
A genomic encoder according to claim 53. The reference conversion unit (2019) applies the same reference conversion (300) to all classes and subclasses of data;
The genome encoder according to claim 54. The reference conversion unit (2019) applies different reference conversions (301, 302, 303) to different classes and subclasses of data;
The genome encoder according to claim 54. Further comprising coding means suitable for performing the coding method according to claim 12.
The genome encoder according to claim 54. A genome decoder (218) for decompressing the compressed genome stream (211), wherein the genome decoder (218) comprises:
A demultiplexer (210) for demultiplexing the compressed genomic data and metadata;
Analysis means (212-214) configured to parse the compressed genome stream into genome blocks (215) of descriptors;
One or more block decoders (216-217) configured to decode genomic blocks of the descriptor into sorted reads of a sequence of nucleotides (211);
A genomic data class decoder configured to selectively decode said categorized reads of the sequence of nucleotides on one or more reference sequences to generate an uncompressed read of the sequence of nucleotides; ,including. Further comprising a reference transform decoder (2113) configured to decode the reference transform descriptor (2112) and generate a transformed reference (2114) used by the genomic data class decoder (219).
A genomic decoder according to claim 58. The one or more reference sequences are stored in a compressed genome stream (211);
A genomic decoder according to claim 59. The one or more reference sequences are provided to the decoder via an out of band mechanism;
A genomic decoder according to claim 59. The one or more reference sequences are constructed at a decoder;
A genomic decoder according to claim 59. One or more reference sequences are converted at the decoder by a reference conversion decoder (2113);
A genomic decoder according to claim 59.   A computer readable medium comprising instructions for causing at least one processor to execute the coding method of claim 12.   60. A computer readable medium containing instructions for causing at least one processor to execute the decoding method of claim 59. Support data for storing a genome encoded according to the method of claim 12.

Images