Gene analysis apparatus, gene analysis method, and gene analysis system

Abstract

A gene analysis apparatus according to the present disclosure comprises: a storage device that stores time-series data received from an electrophoresis device that electrophoreses a fluorescently labeled nucleic acid to acquire the time-series data represented as waveforms of the fluorescence intensity of the nucleic acid; and a processing device that processes the time-series data to determine the nucleotide sequence of the nucleic acid. The processing device executes: processing for acquiring an observed signal by performing mobility correction on the time-series data; processing for estimating a function of a peak width with respect to a migration time from the observed signal; processing for estimating the number and positions of peaks of each waveform in a plurality of small sections of the observed signal; processing for determining peak heights at the respective peak positions on the basis of the area in the vicinity of each peak position of the observed signal; processing for generating a reconstructed signal represented as Gaussian waveforms having the respective peak heights; processing for calculating evaluation values on the basis of the differences between the observed signal and the reconstructed signal; and processing for determining an optimal peak position from among the peak positions on the basis of the evaluation values.

Description

Genetic analysis device, genetic analysis method, and genetic analysis system

[0001] The present disclosure relates to a genetic analysis device, a genetic analysis method, and a genetic analysis system.

[0002] Genetic analyzers (base sequencers) for determining the base sequence of nucleic acids are known. Genetic analyzers are generally called DNA sequencers. Capillary electrophoresis sequencers determine the base sequence by electrophoresing a nucleic acid sample. Capillary electrophoresis sequencers observe the wavelength spectrum of fluorescence emitted from dyes corresponding to four bases (A (adenine), G (guanine), C (cytosine), and T (thymine)) as a time-series input signal to determine the base sequence. This process of determining the base sequence from a time-series input signal is called "base calling."

[0003] Patent Document 1 discloses that "the base sequence determination apparatus includes: (1) a mobility correction unit that outputs a mobility-corrected signal obtained by mobility-correcting a time-series signal of a wavelength spectrum corresponding to each base; (2) a deconvolution unit that executes the following processes: calculating a deconvoluted signal of the mobility-corrected signal for each of a plurality of parameter candidates for a point spread function; calculating a variance of peak intervals for the calculated deconvoluted signal; identifying a parameter of the point spread function using the calculated variance; and outputting the deconvoluted signal corresponding to the point spread function having the identified parameter as an updated deconvoluted signal; (3) a peak extraction unit that extracts a peak waveform from the updated deconvoluted signal and outputs an updated peak-extracted signal; and (4) a sequence identification unit that inputs the updated peak-extracted signal and determines the base sequence" (see Abstract of Patent Document 1).

[0004] International Publication No. 2017 / 130349

[0005] In the technology of Patent Document 1, the smaller the peak width of the point spread function, the higher the reproducibility of the mobility-corrected signal in the deconvoluted signal, and therefore peaks corresponding to point spread functions with narrow peak widths tend to be extracted. This tendency poses a problem in that, in a waveform in which the same base type is repeated, if the peak resolution of the waveform decreases due to the repeated occurrence of the same base type in a region with a long migration time, the number of bases is likely to be determined to be too large. To address this issue, Patent Document 1 calculates the variance of the peak intervals and adds it to the evaluation function. However, it is not easy to adjust the degree to which the variance of the peak intervals is weighted and added to the evaluation function. This may result in an inappropriate weighting depending on the migration environment, resulting in an erroneous determination of the number of bases.

[0006] Therefore, the present disclosure provides a technology for determining base sequences with high accuracy from time-series data showing the results of nucleic acid electrophoresis.

[0007] In order to solve the above problem, the genetic analysis device of the present disclosure is a genetic analysis device capable of communicating with an electrophoresis device that electrophoreses fluorescently labeled nucleic acids to obtain time series data, which is a waveform of the fluorescence intensity of the nucleic acids. The genetic analysis device includes a storage device that stores the time series data received from the electrophoresis device, and a processing device that processes the time series data to determine the base sequence of the nucleic acids. The processing device is characterized by performing the following processes: obtaining an observed signal by performing mobility correction on the time series data; estimating a peak width function of the observed signal with respect to the electrophoresis time; estimating the number of peaks and peak positions for each waveform within a plurality of small sections of the observed signal; determining the peak heights at each of the peak positions based on the areas surrounding the peak positions of the observed signal; generating a reconstructed signal of a Gaussian waveform having each of the peak heights; calculating an evaluation value based on the difference between the observed signal and the reconstructed signal; and determining an optimal peak position from the peak positions based on the evaluation value.

[0008] Furthermore, the genetic analysis method disclosed herein is a genetic analysis method executed by a processing device of a genetic analysis device capable of communicating with an electrophoresis device that electrophoreses fluorescently labeled nucleic acids to acquire time series data, which is a waveform of the fluorescence intensity of the nucleic acids, and is characterized by including the steps of: receiving the time series data from the electrophoresis device by the processing device; acquiring an observation signal by performing mobility correction on the time series data; estimating a peak width function of the observation signal with respect to the electrophoresis time; estimating the number of peaks and peak positions for each waveform within a plurality of small sections of the observation signal; determining peak heights at each of the peak positions based on the areas surrounding the peak positions of the observation signal; generating reconstructed signals of Gaussian waveforms having each of the peak heights; calculating an evaluation value based on the difference between the observation signal and the reconstructed signal; and determining an optimal peak position from among the peak positions based on the evaluation value.

[0009] a processing unit configured to process the time series data to determine the base sequence of the nucleic acid; a processing unit configured to perform mobility correction on the time series data; a processing unit configured to perform a process of estimating a peak width function of the observed signal with respect to electrophoresis time; a processing unit configured to estimate the number of peaks and the peak positions for each of the waveforms within a plurality of small sections of the observed signal; a processing unit configured to generate a reconstructed signal of a Gaussian waveform having each of the peak heights; a processing unit configured to calculate an evaluation value based on the difference between the observed signal and the reconstructed signal; and a processing unit configured to determine an optimal peak position from among the peak positions based on the evaluation value.

[0010] Further features related to the present disclosure will become apparent from the description of this specification and the accompanying drawings. Also, aspects of the present disclosure are achieved and realized by the elements and combinations of various elements and the aspects of the following detailed description and the appended claims. The description of this specification is merely exemplary and does not limit the scope or application of the claims of the present disclosure in any way.

[0011] According to the present disclosure, it is possible to determine the base sequence with high accuracy from time-series data showing the results of electrophoresis of nucleic acids. Problems, configurations, and effects other than those described above will become clear from the following description of the embodiments.

[0012] 1 is a block diagram showing an example of the configuration of a genetic analysis system according to the present disclosure. FIG. 1 is a schematic diagram showing an example of the configuration of an electrophoresis apparatus. FIG. 2 is a flowchart showing an overview of processing performed by a genetic analysis system. FIG. 3 is a flowchart showing electrophoresis processing of a sample by an electrophoresis apparatus. FIG. 4 is a flowchart showing base calling processing by a base calling unit. FIG. 5 is a diagram for explaining an example of block division of time-series data. FIG. 6 is a flowchart showing a method of detecting peaks within a block. FIG. 7 is a diagram for explaining details of processing for generating a reconstructed signal. FIG. 8 is a diagram for explaining details of processing for generating a reconstructed signal. FIG. 9 is a diagram for explaining details of processing for generating a reconstructed signal. FIG. 10 is a diagram for explaining details of processing for generating a reconstructed signal. FIG. 11 is a diagram for explaining the effect of processing for obtaining a reconstructed signal based on an area around a peak position. FIG. 12 is a diagram for explaining the effect of processing for obtaining a reconstructed signal based on an area around a peak position. FIG. 13 is a diagram for explaining a method of correcting the number of peaks. FIG. 14 is a diagram for explaining a method of correcting the number of peaks. FIG. 15 is a diagram for explaining an example of an effect of the present disclosure. FIG. 16 is a diagram for explaining an example of an effect of the present disclosure. FIG. 17 is a diagram for explaining an example of an effect of the present disclosure. FIG. 18 is a diagram for explaining an example of a screen showing the results of base calling. FIG. 19 is a diagram for explaining an example of a GUI screen for correcting parameters related to calculation of the area of ​​an observed signal.

[0013] Hereinafter, embodiments of the present disclosure will be described with reference to the drawings.

[0014] 1 is a block diagram showing an example of the configuration of a genetic analysis system 1 according to a first embodiment. The genetic analysis system 1 includes a data analysis device 100 (genetic analysis device) and an electrophoresis device 200. The data analysis device 100 and the electrophoresis device 200 are connected to each other via a communication cable or wirelessly so that they can communicate with each other.

[0015] The data analysis device 100 includes a central control unit 101, a user interface unit 102, and a memory unit 103. The central control unit 101 controls the operation of the electrophoresis device 200 and performs data processing. The central control unit 101 is, for example, a central processing unit (CPU) and a graphics processing unit (GPU). The memory unit 103 stores programs executed by the central control unit 101, setting information for the electrophoresis device 200, information used for various processes, etc. The memory unit 103 is, for example, a memory and a storage device. The user interface unit 102 is an interface that connects to an input device and an output device, or an interface that connects to an external device via a network. The data analysis device 100 presents information to a user via the user interface unit 102 and also accepts information input by the user.

[0016] The central control unit 101 executes programs stored in the memory unit 103 to operate as a sample information setting unit 104, an electrophoresis device control unit 105, a fluorescence intensity calculation unit 106, and a base call unit 107. In the following description, when processing is described using these functional units as the subject, it indicates that the central control unit 101 is executing the programs.

[0017] The sample information setting unit 104 sets information about the sample (e.g., DNA fragments). The electrophoresis device control unit 105 controls the operation of the electrophoresis device 200 to control the electrophoresis of the sample. The electrophoresis device 200 electrophoreses the sample and acquires electrophoresis data. The electrophoresis data is time-series data of the brightness values ​​(fluorescence intensity) of DNA fragments labeled with fluorescent dyes. The fluorescence intensity calculation unit 106 acquires time-series data indicating the results of the electrophoresis from the electrophoresis device 200. The time-series data includes multiple fluorescence intensity data corresponding to multiple bases.

[0018] The base calling unit 107 analyzes the base sequence of a sample from the time-series data. The base calling unit 107 includes an analysis interval detection unit 108, a mobility correction unit 109, and a base sequencing unit 110. The analysis interval detection unit 108 detects an analysis interval to be analyzed from the time-series data. The mobility correction unit 109 corrects differences in mobility for each base in the time-series data. The base sequencing unit 110 determines the base sequence from the time-series data after mobility correction. Details of the processing by the base sequencing unit 110 will be described later.

[0019] 2 is a schematic diagram showing an example of the configuration of an electrophoresis apparatus 200. The electrophoresis apparatus 200 includes a capillary array 201, a pump mechanism 203, a high-voltage power supply 204, a first ammeter 205, a block 207, a polymer container 209, an anode buffer container 210, an anode electrode 211, a second ammeter 212, a light source 214, an optical detector 215, a detection unit 216, a diffraction grating 217, a thermostatic bath 218, a conveyor 225, a hollow electrode 226, and a load header 229.

[0020] The capillary array 201 includes multiple (e.g., eight) capillaries 202. When a capillary 202 is damaged or its quality deteriorates, it can be replaced with a new capillary array 201. The capillary 202 is composed of a glass tube with an inner diameter of several tens to several hundreds of microns and an outer diameter of several hundreds of microns. The surface of the capillary 202 is coated with a polyimide film to improve its strength. However, in the detection unit 216 where excitation light is irradiated, the polyimide film is removed to allow the emitted light from inside the capillary 202 to easily leak to the outside. The capillary 202 is filled with a separation medium to impart a difference in migration speed during electrophoresis. Separation media are available in both fluid and non-fluid forms, but in the first embodiment, a fluid polymer is used.

[0021] The high-voltage power supply 204 applies a high voltage to the capillary 202. The first ammeter 205 detects the current generated from the high-voltage power supply 204. The second ammeter 212 detects the current flowing through the anode electrode 211.

[0022] The light source 214, optical detector 215, detection unit 216, and diffraction grating 217 constitute an optical detection unit that detects information light obtained from the sample. The light source 214 irradiates the detection unit 216 with excitation light (e.g., laser light). The optical detector 215 detects light emitted from the detection unit 216. The detection unit 216 is a component that acquires sample-dependent information. When detecting a sample in the capillary 202 separated by electrophoresis, the light source 214 irradiates the detection unit 216 with excitation light, thereby generating fluorescence having a wavelength dependent on the sample as information light. Furthermore, the diffraction grating 217 disperses the information light in the wavelength direction, and the optical detector 215 detects the dispersed information light. A user can use the data analysis device 100 to control various functions of the electrophoresis device 200 and acquire electrophoresis data detected by the optical detector 215.

[0023] The hollow electrodes 226 and capillary cathode ends 227 are shown in an enlarged cross-sectional view of the portion circled by a dotted line in Figure 2. The capillary cathode ends 227 are fixed through metal hollow electrodes 226, and the tip of the capillary 202 protrudes from the hollow electrode 226 by about 0.5 mm. All of the hollow electrodes 226 provided on each capillary 202 are mounted integrally on a load header 229. All of the hollow electrodes 226 are electrically connected to a high-voltage power supply 204 mounted in the main body of the apparatus, and function as cathode electrodes when voltage application is required for electrophoresis, sample introduction, etc.

[0024] The end (other end) of the capillary array 201 opposite to the capillary cathode end 227 is bundled together by a capillary head 233. The capillary head 233 can be connected to the block 207 in a pressure-tight, airtight manner. A high voltage from a high-voltage power supply 204 is applied between the load header 229 and the capillary head 233.

[0025] The pump mechanism 203 is composed of a syringe 206 and a mechanism for pressurizing the syringe 206, and injects a polymer into the capillary 202. A block 207 is a connection portion for communicating the syringe 206, the capillary array 201, the anode buffer container 210, and the polymer container 209. The syringe 206 fills new polymer into the capillary 202 from the capillary head 233. The polymer in the capillary 202 is refilled for each measurement to improve measurement performance.

[0026] The thermostatic bath 218 is covered with a heat insulating material to keep the capillaries 202 in the thermostatic bath 218 at a constant temperature, and the temperature is controlled by a heating and cooling mechanism 220. A fan 219 circulates and agitates the air in the thermostatic bath 218, keeping the temperature of the capillary array 201 uniform and constant across its position.

[0027] The transporter 225 transports various containers to the capillary cathode end 227. The transporter 225 is equipped with three electric motors and linear actuators, and is movable in three axial directions: up and down, left and right, and front and back. At least one container can be placed on the moving stage 230 of the transporter 225. The moving stage 230 is also equipped with an electric grip 231, which can grip and release each container. Therefore, the buffer container 221, the washing container 222, the waste container 223, and the sample plate 224 can be transported to the capillary cathode end 227 as needed. Unnecessary containers are stored in a designated storage location within the electrophoresis apparatus 200.

[0028] The electrophoresis device 200 may include sensors for acquiring information about the observation environment that affects electrophoresis (observation environment information). The electrophoresis device 200 in FIG. 2 includes an in-device sensor 240, a polymer sensor 241, and a buffer solution sensor 242.

[0029] The internal sensor 240 is a sensor for acquiring information about the internal environment of the electrophoresis device 200. The internal sensor 240 is, for example, a temperature sensor, a humidity sensor, and an air pressure sensor inside the electrophoresis device 200.

[0030] The polymer sensor 241 is a sensor for acquiring information about the quality of the polymer. The polymer sensor 241 is, for example, a pH sensor, an electrical conductivity sensor, etc. In FIG. 2 , the polymer sensor 241 is installed inside the polymer container 209, but the installation location is not limited thereto.

[0031] The buffer solution sensor 242 is a sensor for obtaining information about the quality of the buffer solution. The buffer solution sensor 242 is, for example, a temperature sensor. Although the buffer solution sensor 242 is installed in the anode buffer container 210 in FIG. 2 , the installation location is not limited thereto. The buffer solution sensor 242 may be installed in the buffer container 221, for example.

[0032] <Outline of Genetic Analysis Flow> FIG. 3 is a flowchart showing an outline of the processing executed by the genetic analysis system 1.

[0033] Step S301: The electrophoresis apparatus control unit 105 of the data analysis apparatus 100 transmits an operation instruction to the electrophoresis apparatus 200. The electrophoresis apparatus 200 then performs an electrophoresis process on the sample to be analyzed. The details of the electrophoresis process will be described later.

[0034] - Step S302 The fluorescence intensity calculation unit 106 of the data analysis device 100 performs spectrum correction to correct the wavelength characteristics of the device.

[0035] Step S303: The fluorescence intensity calculation unit 106 of the data analysis device 100 executes a fluorescence intensity calculation process using the electrophoresis data. Specifically, the fluorescence intensity calculation unit 106 calculates time-series data of the fluorescence intensity of the fluorescent dye from the electrophoresis data, and detects the center position, height, width, etc. of the peak from the time-series data of the fluorescence intensity.

[0036] - Step S304 The mobility correction unit 109 of the data analysis device 100 executes mobility correction processing on the time series data of the fluorescence intensity.

[0037] - Step S305 The base calling unit 107 of the data analysis device 100 executes base calling using the time-series data of the fluorescence intensity corrected based on the result of the mobility correction process, and identifies the base sequence of the sample.

[0038] 4 is a flowchart showing the sample electrophoresis process in step S301 by the electrophoresis device 200. The basic steps of electrophoresis can be broadly divided into sample preparation (step S401), analysis start event (step S402), loading of the migration medium (step S403), preliminary migration (step S404), sample introduction (step S405), migration analysis (step S406), and migration analysis end (step S407).

[0039] Step S401: The user of the electrophoresis device 200 loads samples and reagents into the device as sample preparation before starting analysis. More specifically, first, the buffer container 221 and the anode buffer container 210 are filled with a buffer solution that forms part of the current path. The buffer solution is, for example, an electrolyte solution commercially available from various companies for electrophoresis. The sample to be analyzed is then dispensed into the wells of the sample plate 224. The sample is, for example, a DNA PCR product. A cleaning solution for cleaning the capillary cathode end 227 is then dispensed into the cleaning container 222. The cleaning solution is, for example, pure water. The syringe 206 is then filled with a migration medium for electrophoresis of the sample. The migration medium is, for example, a polyacrylamide-based separation gel or polymer commercially available from various companies for electrophoresis. Furthermore, if deterioration of the capillaries 202 is expected or if the length of the capillaries 202 is to be changed, the capillary array 201 is replaced.

[0040] At this time, the samples loaded onto the sample plate 224 include the actual DNA sample to be analyzed, as well as a positive control, a negative control, and an allelic ladder, each of which is electrophoresed in a different capillary 202. The positive control is, for example, a PCR product containing known DNA and is a control sample used to verify that the DNA is properly amplified by PCR. The negative control is a PCR product containing no DNA and is a control sample used to verify that the PCR amplification product is not contaminated with the user's DNA or dust. The allelic ladder is an artificial sample containing many alleles that may commonly be contained in DNA markers and is typically provided by reagent manufacturers as part of a reagent kit for DNA identification. The allelic ladder is used to fine-tune the correspondence between the DNA fragment length of each DNA marker and the allele. In addition, known DNA fragments labeled with specific fluorescent dyes, called size standards, are mixed with all of the above-mentioned actual sample, positive control, negative control, and allelic ladder samples. The type of fluorescent dye assigned to the size standard varies depending on the reagent kit used.

[0041] The user sets the type of allelic ladder, the type of size standard, the type of fluorescent reagent, and the type of sample set in the wells on the sample plate 224 corresponding to each capillary. In this embodiment, the sample type is specified as any of the following: real sample, positive control, negative control, and allelic ladder. This information is set in the sample information setting unit 104 on the data analysis device 100 via the user interface unit 102.

[0042] Step S402: After completing the sample preparation as described above, the user operates the user interface unit 102 on the data analysis device 100 to instruct the start of analysis. This instruction to start analysis is passed to the electrophoresis device control unit 105. The electrophoresis device control unit 105 then sends an analysis start signal to the electrophoresis device 200, thereby starting the analysis.

[0043] Step S403: The electrophoresis device 200 starts loading the migration medium. This step may be performed automatically after the start of analysis, or may be performed sequentially in response to a control signal sent from the electrophoresis device control unit 105. Loading the migration medium is a procedure for filling the capillaries 202 with new migration medium to form a migration path.

[0044] In the loading of the electrophoretic medium in this embodiment, first, the waste liquid container 223 is carried by the conveyor 225 to a position directly below the load header 229, and the electromagnetic valve 213 is closed so that the used electrophoretic medium discharged from the capillary cathode end 227 can be received. Then, the syringe 206 is driven to load the capillary 202 with new electrophoretic medium, and the used electrophoretic medium is discarded in the waste liquid container 223. Finally, the capillary cathode end 227 is immersed in a cleaning solution in a cleaning container 222 to clean the capillary cathode end 227 contaminated with the electrophoretic medium.

[0045] Step S404: The electrophoresis device 200 performs a preliminary run. This step may be performed automatically or sequentially in response to a control signal sent from the electrophoresis device control unit 105. The preliminary run is a procedure in which a predetermined voltage is applied to the migration medium to prepare the migration medium for electrophoresis. In this embodiment, the preliminary run first involves immersing the capillary cathode end 227 in the buffer solution in the buffer container 221 by the transport device 225 to form a current path. Then, the high-voltage power supply 204 applies a voltage of several kV to several tens of kV to the migration medium for several minutes to several tens of minutes to prepare the migration medium for electrophoresis. Finally, the capillary cathode end 227 is immersed in a cleaning solution in the cleaning container 222 to clean the capillary cathode end 227 contaminated by the buffer solution.

[0046] Step S405: The electrophoresis apparatus 200 introduces a sample. This step may be performed automatically or sequentially in response to a control signal sent from the electrophoresis apparatus control unit 105. In the sample introduction step, sample components are introduced into the migration path. In this embodiment, the transporter 225 first immerses the capillary cathode end 227 in the sample held in the well of the sample plate 224, and then the solenoid valve 213 is opened. This forms a current path, enabling the introduction of sample components into the migration path. Then, the high-voltage power supply 204 applies a pulse voltage to the current path, thereby introducing the sample components into the migration path. Finally, the capillary cathode end 227 is immersed in a cleaning solution in the cleaning container 222 to clean the capillary cathode end 227 contaminated by the sample.

[0047] Step S406: The electrophoresis device 200 performs electrophoretic analysis of the sample. This step may be performed automatically or sequentially in response to control signals sent from the electrophoresis device control unit 105. In electrophoretic analysis, each sample component contained in the sample is separated and analyzed by electrophoresis. In the electrophoretic analysis of this embodiment, first, the capillary cathode end 227 is immersed in the buffer solution in the buffer container 221 by the conveyor 225 to form a current path. Next, the high-voltage power supply 204 applies a high voltage of approximately 15 kV to the current path, generating an electric field in the electrophoresis path. The generated electric field causes each sample component in the electrophoresis path to move toward the detection unit 216 at a speed dependent on the properties of each sample component. In other words, the sample components are separated based on the difference in their migration speeds. The sample components are then detected in order, beginning with the sample components that reach the detection unit 216. For example, if a sample contains many DNA fragments with different base lengths, differences in migration speed will occur depending on the base length, and DNA fragments will arrive at the detection unit 216 in order of shortest base length. Each DNA fragment is labeled with a fluorescent dye that depends on its terminal base sequence. When excitation light from the light source 214 is irradiated onto the detection unit 216, information light, i.e., fluorescence having a wavelength dependent on the sample, is generated from the sample and emitted to the outside. This information light is detected by the optical detector 215. During electrophoretic analysis, the optical detector 215 detects this information light at regular time intervals and transmits image data to the data analysis device 100. Alternatively, to reduce the amount of information transmitted, the luminance of only a portion of the image data may be transmitted instead of the image data. For example, luminance values ​​sampled at only wavelength positions at regular intervals may be transmitted for each capillary. This luminance value data represents the spectral waveform of each capillary. This spectral waveform is stored in the memory unit 103.

[0048] - Step S407 Finally, when the electrophoresis device 200 has acquired the planned image data, it stops applying voltage and ends the electrophoresis analysis.

[0049] <Base Call Flow> FIG. 5 is a flowchart showing the base call processing by the base call unit 107 in step S305.

[0050] Step S501 The analysis interval detection unit 108 detects an analysis interval to be analyzed from the time series data of the fluorescence intensity that has been mobility corrected in step S304.

[0051] Step S502: The base sequence determination unit 110 performs initial peak detection for the detected analysis section. The initial peak detection may be simple maximum value detection based on multiple points. A signal intensity threshold may be set for the maximum value. Furthermore, noise may be removed in advance using a known low-pass filter.

[0052] Step S503: The base sequence determination unit 110 estimates a peak width function using the peak detection results obtained in step S502. The peak width function is a function that represents the relationship between electrophoresis time and peak width. Specifically, the base sequence determination unit 110 first extracts single-base locations where the same base type is not consecutive from the peak detection results, and calculates the peak width from the peak waveform at the single-base location. A known fitting process using a Gaussian function can be used to calculate the peak width. In this embodiment, the standard deviation σ of the fitted Gaussian function is used as the peak width, but the half-width or the like may also be used as the peak width. The base sequence determination unit 110 plots the numerous peak widths thus obtained in a two-dimensional space with the electrophoresis time on the horizontal axis and the peak width on the vertical axis, and obtains a peak width function by polynomial approximation. The above process may be performed for each base type to obtain a peak width function for each base type, or the peak widths of all base types may be plotted in a single space to obtain a peak width function common to all bases. Furthermore, if there is an outlier whose difference from the obtained approximate curve is greater than a predetermined threshold, the outlier may be excluded and the approximate curve may be calculated again. This process may be repeated until there are no more outliers.

[0053] Step S504: The base sequence determination unit 110 estimates a peak interval function using the peak detection results obtained in step S502. The peak interval function is a function that represents the relationship between electrophoresis time and peak intervals. Specifically, the base sequence determination unit 110 calculates the interval between adjacent peaks from the peak detection results. The base sequence determination unit 110 then plots the result in a two-dimensional space with the electrophoresis time on the horizontal axis and the calculated peak interval on the vertical axis, and obtains the peak interval function by polynomial approximation. Furthermore, if there is an outlier whose difference from the obtained approximation curve is greater than a predetermined threshold, the outlier may be excluded, and the approximation curve may be calculated again. This process may be repeated until there are no more outliers.

[0054] Step S505: The base sequence determination unit 110 divides the time-series data into blocks, and performs peak detection (described later) on a block-by-block basis.

[0055] 6 is a diagram illustrating an example of block division of time-series data. As an example, when time-series data is divided by a predetermined threshold value 801, divided sections 802, 803, and 804 are defined as blocks. The threshold value 801 may be determined, for example, as a predetermined percentage of the overall average signal intensity. Alternatively, after locally calculating the average signal intensity, the threshold value 801 may be changed for each migration time by a predetermined percentage for each migration time. This block division is performed for each base type.

[0056] - Step S506 Return to the explanation of Figure 5. In step S506, the base sequence determination unit 110 performs peak detection of the time-series data within each block obtained in step S505. The base sequence determination unit 110 then assigns base sequences to the obtained peak detection positions. Details of step S506 will be described later.

[0057] Step S507: The base sequence determination unit 110 determines whether the processing of step S506 has been performed on all blocks obtained in step S505. If the processing of step S506 has been completed on all blocks (Yes), the process proceeds to step S508. If the processing of step S506 has not been completed on all blocks (No), the process returns to step S506.

[0058] - Step S508: The base sequence determination unit 110 performs final adjustments to the base sequence. Whereas the processing in step S506 is performed for each block of each base type, in step S508, adjustments are performed for the entire base sequence given for all blocks of all base types. The adjustments involve deleting or adding bases based on a threshold value determined in advance for the spacing between adjacent base sequences. This threshold value may be a fixed value, or may be calculated for each electrophoresis time in accordance with the peak width function. Alternatively, bases may be added or deleted after analyzing the peak shape at the corresponding peak position and determining whether or not there is a maximum value. This final adjustment process has the effect of correcting errors in the number of consecutive bases, which are prone to occur when performing step S506 on a block-by-block basis as described above.

[0059] 7 is a flowchart showing the method for detecting peaks within a block in step S506. The peak detection within a block in step S506 includes the following steps S601 to S610. In the following description, the time series data after mobility correction will be referred to as an observed signal in order to distinguish between signals obtained from the time series data at each step.

[0060] Step S601: The base sequence determination unit 110 determines a peak width range. Specifically, the base sequence determination unit 110 estimates a peak width from the migration time of the block position using the peak width function described above, and determines a peak width range that allows a certain fluctuation range from this estimated peak width. An example of the fluctuation range is ±10% of the estimated peak width, but it does not necessarily have to be this value.

[0061] - Step S602 The base sequence determination unit 110 updates the peak width by selecting one peak width from a plurality of peak width candidates at a fixed interval within the range of the peak width determined in step S601.

[0062] Step S603: The base sequence determination unit 110 performs deconvolution processing on the observed signal using the Gaussian function of the peak width updated in step S602 as the point spread function. In this embodiment, the observed signal is modeled by convolution of the original signal with the point spread function. The signal obtained by the deconvolution processing is referred to as the original signal.

[0063] Step S604: The base sequence determination unit 110 estimates peak positions from the original signal obtained by the deconvolution process. One example of a method for estimating peak positions of the original signal is to detect local maxima of the original signal. However, even without using the original signal, one or more candidates for the number of peaks may be determined based on the peak positions before and after the current block and the peak interval estimation function described above, and peak positions may be determined at equal intervals. The base sequence determination unit 110 stores the estimated peak positions in the memory unit 103.

[0064] Step S605: The base sequence determination unit 110 calculates the area of ​​the observed signal around each peak position estimated in step S604. The calculation of the area of ​​the observed signal will be described in detail later. The base sequence determination unit 110 stores the calculated area of ​​the observed signal in the storage unit 103.

[0065] Step S606: The base sequence determination unit 110 generates peak signals that exist only at each peak position based on the area around each peak position obtained in step S605. Here, the peak signal is information that includes the peak position and height. The method for determining the height of the peak signal will be described later. The base sequence determination unit 110 stores the height of the peak signal in the storage unit 103.

[0066] Step S607: The base sequence determination unit 110 generates a reconstructed signal by convolution processing of the peak signal obtained in step S606 with a point spread function. Here, the point spread function is a Gaussian function of the peak width updated in step S602, and is normalized so that its area is 1. Note that in this embodiment, convolution processing of the peak signal obtained in step S606 with the point spread function is used to obtain the reconstructed signal, but it may also be generated by summing up multiple Gaussian functions whose heights are determined from the area.

[0067] Step S608: The base sequence determination unit 110 calculates the difference between the reconstructed signal obtained in step S607 and the observed signal, and sets this as an evaluation value. A typical example of a method for calculating the difference is to calculate the absolute value of the difference between the value of the reconstructed signal and the value of the observed signal at each time, and then calculate the sum of the absolute values ​​within the small interval. The base sequence determination unit 110 stores the evaluation value for the peak width in the storage unit 103.

[0068] Step S609: The base sequence determination unit 110 determines whether the processes from step S602 to step S608 have been performed for all peak width candidates (updated in step S602). If the processes for all peak widths have been completed (Yes), the process proceeds to step S610. If the processes for all peak widths have not been completed (No), the process returns to step S602, and steps S602 to S608 are performed for peak widths for which evaluation values ​​have not been calculated.

[0069] Step S610: The base sequence determination unit 110 selects the peak position corresponding to the peak width with the smallest evaluation value from the evaluation values ​​calculated for all peak widths in step S608. This peak position becomes the base sequence position for the observed signal in each block, as determined in step S506.

[0070] 8A to 8D are diagrams illustrating the details of the reconstructed signal generation process from step S605 to step S607 described above. FIG. 8A illustrates a case where peaks are detected at peak positions P1 and P2 from the observed signal 701 of a certain block in step S604. In this case, the area of ​​the observed signal within a fixed interval centered on peak positions P1 and P2 is calculated. Here, the peak width W at the block position is used as an example of the value of the fixed interval. In FIG. 8A, the peak width W at peak positions P1 and P2 is the same value. Alternatively, separate peak widths W1 and W2 may be used for peak positions P1 and P2. The peak width W is calculated using an estimated peak width function. In FIG. 8A, the area of ​​the observed signal within a range of ±W / 2 (from position B1 to position B2) around peak position P1 is designated A1. Furthermore, the area of ​​the observed signal within a range of ±W / 2 (from position B3 to position B4) around peak position P2 is designated A2.

[0071] In step S607, a peak signal having a height A1 is generated at peak position P1, and a peak signal having a height A2 is generated at peak position P2. The peak signals thus obtained are convolved with a Gaussian function having the peak width set in step S602. Fig. 8B shows an example in which a convolved signal 703 corresponding to the peak signal at peak position P1 and a convolved signal 704 corresponding to the peak signal at peak position P2 are generated, and a reconstructed signal 702 is generated by the sum of these signals.

[0072] FIG. 8C shows an example in which peak positions P1 and P2 are close to each other. In FIG. 8C, because the distance D between peak positions P1 and P2 is D<W, the area range of peak position P1 (positions B1 to B2) and the area range of peak position P2 (positions B3 to B4) partially overlap. In such a case, the area may be calculated while allowing the overlap, or, as shown in FIG. 8D, the area range of peak position P1 and the area range of peak position P2 may be modified using the midpoint B2 (equal to B3) between peak positions P1 and P2 as the boundary. Also, in FIGS. 8A and 8C, positions B1 and B4 are located in areas within the block where the observed signal is non-zero. Alternatively, as shown in FIG. 8D, positions B1 and B4 may be shifted by R1 and R2, respectively, to include the entire observed signal, expanding the area range to a position where the observed signal is zero. This results in areas A1 and A2 becoming closer to the area of ​​the observed signal.

[0073] 9A and 9B are diagrams illustrating the effect of a process for obtaining a reconstructed signal based on the area around a peak position. FIG. 9A shows an example of generating a reconstructed signal 902 using a conventional method when two peak positions are erroneously estimated for an observed signal 901 with one peak. In the conventional method, there is a degree of freedom in the height of the peaks at the peak positions, so even if there are an excessive number of peaks, as shown on the right side of FIG. 9A , the height of any one of the peaks is allowed to be extremely small. This reduces the difference between the reconstructed signal 902 and the observed signal 901, potentially resulting in an erroneous determination that the number of peaks is two.

[0074] 9B shows an example of generating a reconstructed signal 9023 using the method of the present disclosure when two peak positions are erroneously estimated for an observed signal 901 with one peak. By determining the peak height based on the area around the peak position, as in the present disclosure, the height of each peak is constrained as shown on the right side of FIG. 9B, making it less likely that the peak height will be extremely small as in FIG. 9A. This is expected to increase the difference between the reconstructed signal 903 and the observed signal 901, making it more difficult to determine that the number of peaks is two. In other words, by appropriately constraining the peak height based on the area of ​​the observed signal, it is expected that the observed signal will be more easily resolved into an appropriate number of peaks.

[0075] <Correction of the Number of Peaks> In the example described above, the peak position is estimated for each of the multiple peak width candidates in step S604, a reconstructed signal is generated based on the area around the estimated peak position, and the optimal peak position is determined. However, if the number of peaks is estimated incorrectly for all of the peak width candidates, the correct number of peaks cannot be estimated. In such cases, it is desirable to be able to correctly correct the number of peaks.

[0076] 10A to 10D are diagrams illustrating a technique for correcting the number of peaks based on the difference between the reconstructed signal and the observed signal. FIG. 10A illustrates a case where the difference value D at the position where the difference between the observed signal 1001 and the reconstructed signal 1002 is maximum (where the observed signal 1001 is greater than the reconstructed signal 1002) exceeds a predetermined threshold. In this case, as shown in FIG. 10B , a new peak is inserted at the position (P4) of the difference value D. The height of the peak at the peak position P4 can be determined based on the area around the peak position, as described above. The evaluation value described above, obtained from the difference between the reconstructed signal 1003 obtained in this way and the reconstructed signal 1002 without the peak insertion, and the smaller evaluation value can be used.

[0077] FIG. 10C illustrates a case where the difference value D at the position where the difference between the observed signal 1004 and the reconstructed signal 1005 is maximum (where the reconstructed signal 1005 is greater than the observed signal 1004) exceeds a predetermined threshold. In this case, as shown in FIG. 10D , peak P4 located near the position of the difference value D is deleted. The peak closest to the position of the difference value D may be selected, or the smaller of the nearby peaks may be selected. Alternatively, peaks near the position of the difference value D may be deleted one by one, a reconstructed signal 1006 may be generated for each peak, an evaluation value obtained from the difference with the observed signal 1004 may be calculated, and the peak with the smaller evaluation value may be deleted. The evaluation values ​​described above may be calculated for the reconstructed signal 1006 thus obtained and the reconstructed signal 1005 without the peaks deleted, and the smaller evaluation value may be used.

[0078] 11A to 11C are diagrams illustrating an example of the effect of the present disclosure. FIG. 11A shows a case in which an observed signal 1101 of T and an observed signal 1102 of C in the base sequence "TTC" are present, and two noise-induced peaks are observed at the position of the observed signal 1101 corresponding to the second T. When the spacing between such noise-induced peaks is clearly smaller than the spacing between the surrounding peaks, processing can be applied such as treating them as a single peak and excluding one of the peaks. However, when there is a certain distance between the two peaks as in FIG. 11A, conventional technology is unlikely to be able to determine that they are one peak, and will likely treat them as two peaks and output the base sequence "TTTC."

[0079] According to an embodiment of the present disclosure, as shown in FIG. 11B , the height of a peak is determined based on the area surrounding the peak. Therefore, the difference between the reconstructed signal 1103 obtained from nearby peaks and the observed signal 1101 is large, resulting in a high evaluation value and making it difficult to determine the positions of such two peaks. Furthermore, as shown in FIG. 11C , if the maximum difference exceeds a threshold, one of the peaks can be deleted as described above to determine a single peak, and the correct base sequence "TTC" can be output. Thus, according to the technology of the present disclosure, it is possible to correctly output base call results as fewer peaks than expected from the observed signal. Therefore, according to the technology of the present disclosure, peak detection that is robust against noise can be performed.

[0080] <Example of user interface screen> By displaying on the user interface unit 102 how the peak width, peak position, etc. were determined before the individual base sequences (base call results) were output, the user can check the details of the base call results.

[0081] FIG. 12 is a diagram showing an example of a GUI screen 1200 showing the results of base calling. The GUI screen 1200 is displayed on the user interface unit 102 (display device). The GUI screen 1200 has a base call result display section 1201 showing a chromatogram (observed signal) and the finally determined base sequence. When the user clicks on a base that they wish to confirm in the base call result display section 1201, the finally determined peak width, peak position, reconstructed signal, and evaluation value for that block are presented in display section 1202. The solid line in the chromatogram in display section 1202 indicates the observed signal, and the dotted line indicates the reconstructed signal. As in display section 1203, the peak position, reconstructed signal, and evaluation value may be displayed for each of all peak width candidates (σ1 to σ6). This allows the user to confirm whether the adopted peak width candidate is appropriate.

[0082] The GUI screen 1200 shown in FIG. 12 can be used, for example, for operations checks by a user or service person of the gene analysis system 1.

[0083] FIG. 13 shows an example of a GUI screen 1300 for modifying parameters related to area calculation of observed signals and performing base calling again. The GUI screen 1300 includes the same base call result display section 1201 as in FIG. 12 and a modification operation section 1301. When a user clicks on a base that the user wants to confirm in the base call result display section 1201, the modification operation section 1301 displays the estimated peak position in that block, the width of the Gaussian waveform for area calculation (peak width), the height of the peak signal, the number of peaks, and the like. The modification operation section 1301 is configured to allow the user to modify this information. The modification operation section 1301 in FIG. 13 shows the user moving the width for area calculation (boundary of the area range) by operating the pointer. When the user modifies the width for area calculation, a peak with a height corresponding to the modified value is calculated.

[0084] Summary of First Embodiment As described above, the genetic analysis system 1 according to the first embodiment includes an electrophoresis apparatus 200 that electrophoreses fluorescently labeled nucleic acids to acquire time-series data representing the fluorescence intensity waveforms of the nucleic acids, and a data analysis apparatus 100 that can communicate with the electrophoresis apparatus 200. The central control unit 101 of the data analysis apparatus 100 executes the following processes: acquiring an observed signal by performing mobility correction on the time-series data; estimating a peak width function for the observed signal relative to the electrophoresis time; estimating the number and position of peaks for each waveform within multiple blocks (subsections) of the observed signal; determining peak heights at each peak position based on the area surrounding the peak position of the observed signal; generating a reconstructed signal of a Gaussian waveform having each peak height; calculating an evaluation value based on the difference between the observed signal and the reconstructed signal; and determining an optimal peak position from the evaluation value. This allows the number of consecutive bases to be detected with high accuracy from the time-series data representing the results of electrophoresis, thereby enabling highly accurate base sequence determination.

[0085] [Modifications] The present disclosure is not limited to the above-described embodiments and includes various modifications. For example, the above-described embodiments have been described in detail to clearly explain the present disclosure, and it is not necessary to include all of the described configurations. Furthermore, a part of one embodiment can be replaced with a configuration of another embodiment. Furthermore, a configuration of another embodiment can be added to a configuration of one embodiment. Furthermore, a part of the configuration of each embodiment can be added to, deleted from, or substituted for a part of the configuration of another embodiment.

[0086] DESCRIPTION OF SYMBOLS 1... Genetic analysis system 100... Data analysis device (genetic analysis device) 200... Electrophoresis device 101... Central control unit (processing device) 102... User interface unit (display device, input device) 103... Memory unit (memory device) 104... Sample information setting unit 105... Electrophoresis device control unit 106... Fluorescence intensity calculation unit 107... Base calling unit 108... Analysis interval detection unit 109... Mobility correction unit 110... Base sequence determination unit

Claims

1. A genetic analysis device capable of communicating with an electrophoresis device that electrophoreses fluorescently labeled nucleic acids to obtain time-series data that is a waveform of the fluorescence intensity of the nucleic acids, the genetic analysis device comprising: a storage device that stores the time-series data received from the electrophoresis device; and a processing device that processes the time-series data to determine the base sequence of the nucleic acids, the processing device performing the following processes: obtaining an observed signal by performing mobility correction on the time-series data; estimating a peak width function of the observed signal with respect to electrophoresis time; estimating the number of peaks and peak positions for each waveform in a plurality of small sections of the observed signal; determining the peak heights at each of the peak positions based on the areas surrounding the peak positions of the observed signal; generating a reconstructed signal of a Gaussian waveform having each of the peak heights; calculating an evaluation value based on the difference between the observed signal and the reconstructed signal; and determining an optimal peak position from among the peak positions based on the evaluation value.

2. The genetic analysis device described in claim 1, characterized in that, in the process of determining the peak height at each of the peak positions, the processing device obtains the peak interval at the peak positions based on a peak interval function estimated in advance, and determines the peak height at each of the peak positions based on the area of ​​the observed signal in the section of the peak interval centered at the peak position.

3. The genetic analysis device described in claim 2, characterized in that the processing device extends the section for obtaining the area to include all of the observed signals, and determines the peak height based on the area of ​​the observed signals in the extended section.

4. The genetic analysis device of claim 2, wherein the processing device, when there are two or more peak positions within the small section, divides the area of ​​the observed signal between the two or more peak positions, and determines the peak height at each of the peak positions based on the divided area.

5. The genetic analysis device of claim 1, characterized in that the processing device further performs at least one of the following processes: a process of deleting peaks at locations where the value of the reconstructed signal is greater than a predetermined threshold value relative to the value of the observed signal; and a process of inserting peaks at locations where the value of the reconstructed signal is smaller than a predetermined threshold value relative to the value of the observed signal, and further performs a process of calculating the evaluation value again.

6. The genetic analysis device of claim 1, characterized in that, in the process of estimating the number of peaks and the peak positions for each waveform within a plurality of small sections of the observed signal, when the observed signal includes a small section having two or more peaks spaced apart at intervals smaller than the peak spacing estimated from the observed signal, the processing device determines a number of peak positions in the small section that is smaller than the number of peaks spaced apart.

7. The genetic analysis device of claim 1, further comprising a display device, wherein the processing device further executes a process of displaying at least one of the base sequence of the nucleic acid, the optimal peak position, the observed signal, the reconstructed signal, and the evaluation value on the display device.

8. The genetic analysis device according to claim 7, characterized in that, in the process of determining the peak height at each of the peak positions, the processing device further executes the process of obtaining the peak interval at the peak positions based on a peak interval function estimated in advance, determining the peak height at each of the peak positions based on the area of ​​the observed signal in the section of the peak interval centered at the peak position, and displaying at least one of the peak interval, the peak positions, and the number of peaks on the display device so that the user can modify it.

9. A genetic analysis method executed by a processing device of a genetic analysis device capable of communicating with an electrophoresis device that electrophoreses fluorescently labeled nucleic acids and acquires time-series data that is a waveform of the fluorescence intensity of the nucleic acids, the genetic analysis method comprising the steps of: receiving the time-series data from the electrophoresis device; acquiring an observation signal by performing mobility correction on the time-series data; estimating a peak width function of the observation signal with respect to the electrophoresis time; estimating the number of peaks and peak positions for each waveform in a plurality of small sections of the observation signal; determining peak heights at each of the peak positions based on the areas surrounding the peak positions of the observation signal; generating reconstructed signals of Gaussian waveforms having each of the peak heights; calculating an evaluation value based on the difference between the observation signal and the reconstructed signal; and determining an optimal peak position from among the peak positions based on the evaluation value.

10. The genetic analysis method of claim 9, characterized in that determining the peak height at each of the peak positions includes: obtaining the peak interval at the peak positions based on a peak interval function estimated in advance; and determining the peak height at each of the peak positions based on the area of ​​the observed signal in the section of the peak interval centered at the peak position.

11. The genetic analysis method described in claim 10, characterized in that determining the peak height at each of the peak positions includes: expanding the section for obtaining the area to include all of the observed signals; and determining the peak height based on the area of ​​the observed signals in the expanded section.

12. The genetic analysis method described in claim 10, characterized in that determining the peak height at each of the peak positions includes, when there are two or more peak positions within the small section, dividing the area of ​​the observed signal between the two or more peak positions, and determining the peak height at each of the peak positions based on the divided area.

13. The genetic analysis method of claim 9, further comprising the processing device at least one of: deleting peaks at locations where the value of the reconstructed signal is greater than a predetermined threshold value relative to the value of the observed signal; and inserting peaks at locations where the value of the reconstructed signal is smaller than a predetermined threshold value relative to the value of the observed signal; and recalculating the evaluation value.

14. The genetic analysis method of claim 9, characterized in that estimating the number of peaks and the peak positions for each waveform within a plurality of small sections of the observed signal includes, when the observed signal includes a small section having two or more peaks spaced apart at intervals smaller than the peak spacing estimated from the observed signal, determining a smaller number of peak positions in the small section.

15. The genetic analysis method of claim 9, further comprising the steps of: the genetic analysis device further comprising a display device; and the processing device further comprising causing the display device to display at least one of the base sequence of the nucleic acid, the optimal peak position, the observed signal, the reconstructed signal, and the evaluation value.

16. The genetic analysis method of claim 15, wherein determining the peak height at each of the peak positions includes: obtaining the peak interval at the peak positions based on a peak interval function estimated in advance; and determining the peak height at each of the peak positions based on the area of ​​the observed signal in the section of the peak interval centered at the peak position; and wherein the genetic analysis method further includes displaying at least one of the peak interval, the peak positions, and the number of peaks on the display device so that the user can modify it.

17. A genetic analysis system comprising: an electrophoresis apparatus that electrophoreses fluorescently labeled nucleic acids and detects fluorescence from the nucleic acids to obtain time-series data that is a waveform of fluorescence intensity; and a genetic analysis apparatus configured to be able to communicate with the electrophoresis apparatus, wherein the genetic analysis apparatus has: a storage device that stores the time-series data received from the electrophoresis apparatus; and a processing device that processes the time-series data to determine the base sequence of the nucleic acid, wherein the processing device performs the following processes: obtaining an observed signal by performing mobility correction on the time-series data; estimating a peak width function of the observed signal with respect to electrophoresis time; estimating the number of peaks and peak positions for each waveform in a plurality of small sections of the observed signal; determining the peak heights at each of the peak positions based on the areas surrounding the peak positions of the observed signal; generating a reconstructed signal of a Gaussian waveform having each of the peak heights; calculating an evaluation value based on the difference between the observed signal and the reconstructed signal; and determining an optimal peak position from among the peak positions based on the evaluation value.

Images