JP4805586B2 - Gene structure prediction method and gene structure prediction program - Google Patents

Description

  The present invention relates to a gene structure prediction method and a gene structure prediction program, and in particular, is measured with the tiling array using a tiling array in which partial base sequences regularly extracted from genomic base sequences are arranged as probes. Further, the present invention relates to a gene structure prediction method and a gene structure prediction program for predicting an unknown gene region and structure from the genome base sequence based on tiling array data relating to the expression intensity of each probe.

  So far, for example, Non-Patent Documents 1 to 3 have been disclosed as techniques for predicting the structure of a gene. According to the techniques described in Non-Patent Documents 1 to 3, it was possible to predict the structure of a known gene from only the base sequence data.

  Here, the DNA microarray is very useful as a technique for simultaneously measuring the expression levels of many genes as fluorescence intensity (expression intensity). In particular, tiling arrays are very useful because they can simultaneously measure the expression level of the entire region of the genomic base sequence as fluorescence intensity (expression intensity) (see, for example, Non-Patent Documents 4 to 8). In the tiling array, a plurality of probes are arranged very densely and regularly along the base sequence of the genome regardless of whether or not they are gene regions. As a result, it was possible to discover the presence of an unknown gene using only tiling array data relating to the expression intensity of each probe measured with a tiling array. In the tiling array, a gene region having an average length corresponds to several tens to several hundreds of probes.

Lukashin, AV. Borodovsky, M .; “GeneMark. Hmm: new solutions for gene finding.”, Nucleic Acids Res. , 26, p. 1107-1115, 1998 Uberbacher, EC. , Xu, Y. , Mural, RJ. , “Discovering and understanding genes in human DNA sequence using GRIL”, Methods Enzymol. , 266, p. 259-281, 1996 Snyder, EE. , Stormo, GD. , “Identification of coding regions in genomic DNA sequences: an application of dynamic programming and neural networks.”, Nucleic Acids Res. , 21, p. 607-613, 1993 Shomaker D. D. et al. "Experimental annotation of the human genome using microarray technology.", Nature, 409, p. 922, 2001 Kaplanov P.M. et al. , “Large-scale transcriptional activity in chromosomes 21 and 22.”, Science, 296, p. 916-919, 2002 Rinn, J.M. L. et al. "The transcribable activity of human chromosome 22.", Genes Dev. , 17, p. 529-540, 2003 Yamada K.K. et al. , “Empirical Analysis of Transcriptional Activity in the Arabidopsis Genome.”, Science, 302, p. 842-845, 2003 Stolc V.E. et al. "A gene expression map for the euchromatic genome of Drosophila melanogaster.", Science, 306, p. 655-660, 2004

  However, since tiling array data generally contains a lot of noise, it is not possible to clearly distinguish exons and introns present in the genome base sequence only from tiling array data. Specifically, in the case of a method for determining whether or not an exon is based on whether or not the expression intensity measured by the tiling array exceeds the threshold, even if the threshold is adjusted, it is below the threshold within the exon region. There is an expression intensity, or there is an expression intensity exceeding a threshold outside the exon region. Therefore, there was a problem that even though the presence of an unknown gene could be predicted based on tiling array data, the region and structure of the gene could not be predicted.

  In addition, in the conventional technique for predicting the structure of a gene, even in the prediction process, even a gene expression region that is statistically significant may be determined not to be a gene region. There was a problem that it was not always possible to predict with sufficient accuracy.

  The present invention has been made in view of the above problems, and provides a gene structure prediction method and a gene structure prediction program capable of accurately predicting an unknown gene region and structure from a genomic base sequence based on tiling array data. The purpose is to provide.

  Conventional gene structure prediction methods that use a Markov model are often used (for example, “Asai K., Hayamizu, S. and Handa, K.,“ Prediction of protein secondary structure by the label ”. “Computer Applications for Biosciences, 9, p. 141-146, 1993” and “Krogh, A., Brown, M., Mian, IS. Sjorander, K., Haussler, D.,“ Hidden Markovel. to protein model ng. ”, J. Mol. Biol., 235, p. 1501-1531, 1994”, “Tanaka, H., Ishikawa, M., Asai, K., Konagaya, A.,“ Hidden Markov model and : Study of the equivalence and possibilities. ", International Conference on Intelligent Systems for Molecular Biology, 1, p. 395-401, 1993"). However, in the conventional gene structure prediction method, for example, no Markov model has been applied in consideration of the likelihood calculated from tiling array data. In the first place, no attempt has been made to predict the structure of a gene in consideration of tiling array data.

Therefore, in order to achieve the above object , the gene structure prediction method according to the present invention includes a partial base regularly extracted from a genome base sequence in the gene structure prediction method for predicting the structure of a gene based on the base sequence. Using a tiling array in which the sequences are arranged as probes, and predicting the structure of the gene from the genomic base sequence based on tiling array data on the expression intensity of each probe measured by the tiling array. To do.

The gene structure prediction method according to the present invention uses a tiling array in which partial base sequences regularly extracted from genomic base sequences are arranged as probes, and the expression of each probe measured with the tiling array. Data acquisition step for obtaining tiling array data related to strength and genomic base sequence data related to genomic base sequence, and estimation of gene expression region from genomic base sequence based on tiling array data acquired in the data acquisition step An expression region estimation step and a base corresponding to the expression region estimated in the expression region estimation step as a starting point, for each base that is adjacent to the starting point, a gene based on the tiling array data and the genomic base sequence data By determining attributes related to the structure of the gene, the structure of the gene Characterized in that it comprises a gene structure prediction step of predicting, a.

Moreover, gene structure prediction method according to the present invention, in the gene structure prediction method, the expression region estimation step, the expression determining the expression region candidate is a candidate for the expression region of the gene based on the tiling array data A region candidate determination step, an expression region candidate determination step for determining whether the expression region candidate determined in the expression region candidate determination step is statistically significant, and a determination that the expression region candidate determination step is significant In that case, the method further includes an expression region determination step of determining an expression region candidate as the expression region.

Moreover, gene structure prediction method according to the present invention, in the gene structure prediction method, the expression region candidate determining step, a target region of a predetermined length in the genomic nucleotide sequence, the expression of probes contained in the area A median calculation step for calculating the median value of the intensity; a region moving step for moving the region that is the median calculation target in the median value calculating step along the genome base sequence; the median value calculating step and the region shifting An expression region candidate selection step of selecting a maximum one of a plurality of median values accumulated by repeatedly executing the step, and selecting a region that was a target of calculation of the selected maximum median value as the expression region candidate; , Further included.

Moreover, gene structure prediction method according to the present invention, in the gene structure prediction method, the gene structure prediction step as a target base in which expression region base corresponding to the expression region estimated by the expression region estimating step, An expression region base probability determination step for determining a probability related to the attribute of the expression region base using a Markov model, and an adjacent base adjacent to the expression region base targeted in the expression region probability determination step An adjacent base probability determination step for determining a state transition probability for an attribute using a Markov model, and the adjacent base probability determination step are repeatedly performed up to the end base of the genome base sequence, and the likelihood for the attribute of each base is maximized. A likelihood determining step that uses a likelihood method, and the probability, the state transition probability, and the likelihood An attribute determination step of determining a combination of the attributes of each base in Zui, characterized in that it further comprises a.

Moreover, gene structure prediction method according to the present invention, in the gene structure prediction method, based on the expression intensity of the probe corresponding to exon region included in the structure of the gene predicted by the gene structure prediction step, of the gene An expression value calculating step of calculating an expression value;

Moreover, the gene structure prediction method according to the present invention is the above-described gene structure prediction method, wherein an unpredicted region that is a genomic base sequence excluding the region of the gene whose structure is predicted in the gene structure prediction step is used as a prediction target region. A prediction target region determination step for determining, wherein the expression region estimation step and the gene structure prediction step are performed again on the prediction target region determined in the prediction target region determination step.

Moreover, gene structure prediction method according to the present invention, in the gene structure prediction method, the prediction target region determining step, and the non-prediction region extraction step of extracting the non-prediction region from the genomic nucleotide sequence, wherein the non-prediction region The base length of the unpredicted region extracted in the extraction step is calculated, the base length determination step for determining whether or not the calculated base length is equal to or less than a predetermined base length, and the predetermined base length or less in the base length determination step A non-predicted region extracted in the unpredicted region extraction step, and a prediction target region determining step for determining the unpredicted region as the prediction target region.

The present invention also relates to a gene structure prediction program, and a gene structure prediction program that causes a computer to execute the gene structure prediction method according to the present invention uses a partial base sequence regularly extracted from a genomic base sequence as a probe. Using the arranged tiling array, a data acquisition step for acquiring tiling array data related to the expression intensity of each probe measured with the tiling array and a genomic base sequence data regarding the genomic base sequence, and acquiring in the data acquisition step Based on the tiling array data, the expression region estimation step for estimating the expression region of the gene from the genome base sequence, and the base corresponding to the expression region estimated in the expression region estimation step as a starting point, both adjacent to the starting point For each base connected to By determining the attribute about the structure of the gene on the basis of data and the genome sequence data, characterized in that it contains a gene structure prediction step of predicting a structure of the gene, a.

Moreover, gene structure prediction program according to the present invention, in the gene structure prediction programs, the expression region estimation step, the expression determining the expression region candidate is a candidate for the expression region of the gene based on the tiling array data A region candidate determination step, an expression region candidate determination step for determining whether the expression region candidate determined in the expression region candidate determination step is statistically significant, and a determination that the expression region candidate determination step is significant In that case, the method further includes an expression region determination step of determining an expression region candidate as the expression region.

Moreover, gene structure prediction program according to the present invention, in the gene structure prediction programs, the expression region candidate determining step, a target region of a predetermined length in the genomic nucleotide sequence, the expression of probes contained in the area A median calculation step for calculating the median value of the intensity; a region moving step for moving the region that is the median calculation target in the median value calculating step along the genome base sequence; the median value calculating step and the region shifting An expression region candidate selection step of selecting a maximum one of a plurality of median values accumulated by repeatedly executing the step, and selecting a region that was a target of calculation of the selected maximum median value as the expression region candidate; , Further included.

Moreover, gene structure prediction program according to the present invention, in the gene structure prediction programs, the gene structure prediction step, as the target base at which expression region base corresponding to the expression region estimated by the expression region estimating step, An expression region base probability determination step for determining a probability related to the attribute of the expression region base using a Markov model, and an adjacent base adjacent to the expression region base targeted in the expression region probability determination step An adjacent base probability determination step for determining a state transition probability for an attribute using a Markov model, and the adjacent base probability determination step are repeatedly performed up to the end base of the genome base sequence, and the likelihood for the attribute of each base is maximized. A likelihood determining step determined using a likelihood method, the probability, the state transition probability and the And further comprising a, an attribute determination step of determining a combination of the attributes of each base on the basis of the likelihood.

Moreover, gene structure prediction program according to the present invention, in the gene structure prediction program, based on the expression intensity of the probe corresponding to exon region included in the structure of the gene predicted by the gene structure prediction step, of the gene An expression value calculating step of calculating an expression value;

Moreover, gene structure prediction program according to the present invention, in the gene structure prediction programs, the non-predicted area is the genomic nucleotide sequence structure excluding the region of the gene that is predicted by the gene structure prediction step as a prediction target region A prediction target region determination step for determining, wherein the expression region estimation step and the gene structure prediction step are performed again on the prediction target region determined in the prediction target region determination step.

Moreover, gene structure prediction program according to the present invention, in the gene structure prediction program, the predicted target region determining step includes a non-predictive region extraction step of extracting the non-prediction region from the genomic nucleotide sequence, wherein the non-prediction region The base length of the unpredicted region extracted in the extraction step is calculated, the base length determination step for determining whether or not the calculated base length is equal to or less than a predetermined base length, and the predetermined base length or less in the base length determination step A non-predicted region extracted in the unpredicted region extraction step, and a prediction target region determining step for determining the unpredicted region as the prediction target region.

According to the gene structure prediction method of the present invention, using a tiling array in which partial base sequences regularly extracted from genomic base sequences are arranged as probes, the expression of each probe measured by the tiling array is used. Since the gene structure is predicted from the genome base sequence based on the tiling array data on the strength, the effect that the region and structure of an unknown gene can be accurately predicted from the genome base sequence based on the tiling array data Play.

Further, according to the gene structure prediction method and the gene structure prediction program of the present invention, a tiling array using a partial base sequence regularly extracted from a genomic base sequence as a probe is used. Obtain the tiling array data related to the measured expression intensity of each probe and genomic base sequence data related to the genomic base sequence, and based on the acquired tiling array data, estimate the gene expression region from the genomic base sequence, Predict the structure of a gene by determining the attributes related to the structure of the gene based on the tiling array data and genomic base sequence data for each base that is adjacent to the base corresponding to the estimated expression region. Therefore, based on the tiling array data, An effect that the region and the structure of the gene can be predicted accurately.

In addition, according to the gene structure prediction method and the gene structure prediction program of the present invention, expression region candidates that are gene expression region candidates are determined based on tiling array data, and the determined expression region candidates are statistically determined. If it is determined whether it is significant, and if it is determined to be significant, the expression region candidate is determined as the expression region, so that the expression region can be estimated easily and accurately using existing statistical methods Play.

Further, according to the gene structure prediction method and the gene structure prediction program according to the present invention, for a region of a predetermined length in the genomic base sequence, the median of the expression intensity of the probes contained in the region is calculated, By moving the region for which the value is to be calculated along the genome base sequence and repeatedly executing these processes, the largest one of the accumulated medians is selected and the selected median is calculated. Since the target region is selected as an expression region candidate, there is an effect that an appropriate expression region candidate can be selected in consideration of variations in expression intensity of each probe measured by the tiling array.

Further, according to the gene structure prediction method and the gene structure prediction program according to the present invention, for the expression region base that is the base corresponding to the estimated expression region, the probability relating to the attribute of the expression region base is determined using a Markov model. Determine the state transition probability related to the attribute of the adjacent base using the Markov model for the adjacent base adjacent to the target expression region base, and perform the process of determining the state transition probability at the end of the genome base sequence. After repeatedly executing up to the base, the likelihood regarding the attribute of each base is determined using the maximum likelihood method, and the combination of the attributes of each base is determined based on the determined probability, state transition probability, and likelihood. There is an effect that it is possible to predict the region and structure of the above with higher accuracy.

Further, according to the gene structure prediction method and the gene structure prediction program according to the present invention, the expression value of the gene is calculated based on the expression intensity of the probe corresponding to the exon region included in the predicted gene structure. The effect is that an expression level of an unknown gene can be obtained.

According to the gene structure prediction method and the gene structure prediction program of the present invention, an unpredicted region that is a genomic base sequence excluding a region of a gene whose structure is predicted is determined as a prediction target region, and the determined prediction target Since the process of estimating the expression area and the process of predicting the gene structure are performed again on the region, there is an effect that the regions and structures of a plurality of genes can be efficiently predicted from the genome base sequence.

According to the gene structure prediction method and gene structure prediction program of the present invention, the unpredicted region is extracted from the genome base sequence, the base length of the extracted unpredicted region is calculated, and the calculated base length is a predetermined base. If it is determined whether or not the length is equal to or shorter than the predetermined length, and it is determined that the length is not shorter than a predetermined base length (for example, 1000 base length), the extracted unpredicted region is determined as the prediction target region. As a result, it is possible to eliminate the region, and as a result, it is possible to determine the prediction target region suitable for the structure prediction.

  Embodiments of a gene structure prediction method and a gene structure prediction program according to the present invention will be described below in detail with reference to the drawings. Note that the present invention is not limited to the embodiments.

First, terms used in this specification are defined.
A “gene” is a set of mature mRNA (messenger RNA) or a partial region (exon) of a genomic base sequence corresponding to the mRNA.
The “gene region” is the shortest partial region including all of the set of exons belonging to one gene on the genome base sequence.
“Gene structure” is information indicating exon and intron regions in one gene region.
“Tiling array” refers to a base sequence of a certain length (20 to 100 bases) extracted from a genome base sequence while shifting it at regular intervals (1 to 200 bases), and probes oligonucleotides corresponding to the extracted base sequences. A DNA microarray. The tiling array not only adds an oligonucleotide that perfectly matches the extracted base sequence (perfectly matched oligonucleotide) as a probe, but also an oligonucleotide that does not completely match the extracted base sequence (mismatched oligonucleotide). ) May be added as a comparative probe.
“Expression intensity” is numerical data proportional to the amount of mRNA contained in a sample obtained from a biological tissue, and is measured by the sample hybridizing to a corresponding probe.
“Tiling array data” is data relating to the expression intensity of all probes in the tiling array.
The “significant expression region” means that “the expression intensity of a probe corresponding to a region of a certain base length (for example, 750 base length) in the genome base sequence was measured by chance due to noise”. This is the region where no hypothesis is rejected at a significant level based on statistical tests. The significant expression region is assumed to be a part of the gene region.

  Below, the structure of the gene structure prediction apparatus 100 which is an apparatus for implement | achieving this invention is demonstrated with reference to FIG. FIG. 1 is a block diagram showing a configuration of the gene structure prediction apparatus 100, and conceptually shows only a portion related to the present invention in the configuration.

  As shown in FIG. 1, the gene structure prediction apparatus 100 includes a control unit 102 such as a CPU that comprehensively controls the gene structure prediction apparatus 100, a communication device such as a router, and a wired or wireless communication line such as a dedicated line. A communication interface unit 104 that connects the gene structure prediction apparatus 100 to the network 300 via a network, a storage unit 106 that stores various databases, tables, files, and the like, and an input / output unit that connects to the input device 112 and the output device 114 The interface unit 108 is configured to be communicable via an arbitrary communication path. The network 300 has a function of connecting the gene structure prediction apparatus 100 and the external system 200 so that they can communicate with each other, and is, for example, the Internet or a LAN. Further, the external system 200 is connected to the gene structure prediction apparatus 100 via the network 300 so as to be able to communicate with each other, and provides an external database, genome programs, and the like regarding genome base sequence data, tiling array data, various parameters, and the like. Functions, etc. Each function of the external system 200 is realized by a CPU, a disk device, a memory device, an input device, an output device, a communication control device, and the like in the hardware configuration of the external system 200 and a program for controlling them. The external system 200 may be configured as a WEB server, an ASP server, or the like, and the hardware may be configured by an information processing apparatus such as a commercially available workstation or personal computer, and an accessory device thereof.

  The storage unit 106 is a storage unit. Specifically, a memory device such as a RAM or a ROM, a fixed disk device such as a hard disk, a flexible disk, an optical disk, or the like can be used. The storage unit 106 stores a genomic base sequence data file 106a, a tiling array data file 106b, an expression region data file 106c, a gene structure data file 106d, and an expression value data file 106e as shown in the figure. The genomic base sequence data file 106a stores genomic base sequence data. Specifically, as shown in FIG. 2, a base number for uniquely identifying the position of each base in the genome base sequence and a base symbol (A, T, G, C) representing the base corresponding to the base number Are stored in association with each other. Returning again to FIG. 1, the tiling array data file 106b stores tiling array data. Specifically, as shown in FIG. 3, the probe number for uniquely identifying the probe arranged in the tiling array, the expression intensity of the probe corresponding to the probe number, and the probe range corresponding to the probe number ( A set of a start base number and an end base number in the genome base sequence) are stored in association with each other. Returning to FIG. 1 again, the expression region data file 106c stores data related to the expression region estimated by the expression region estimation unit 102b described later. The gene structure file 106d stores data related to the region and structure of a gene predicted by a gene structure prediction unit 102c described later. The expression value data file 106e stores data relating to expression values calculated by an expression value calculation unit 102d described later.

  The communication interface unit 104 mediates communication between the gene structure prediction device 100 and the network 300 (or a communication device such as a router). That is, the communication interface unit 104 has a function of communicating data with other terminals via a communication line.

  The input / output interface unit 108 is connected to the input device 112 and the output device 114. Here, in addition to a monitor (including a home television), a speaker or a printer can be used as the output device 114 (the output device 114 may be described as a monitor below). In addition to the keyboard, mouse, and microphone, the input device 112 can be a monitor that realizes a pointing device function in cooperation with the mouse.

  The control unit 102 has an internal memory for storing a control program such as an OS (Operating System), a program that defines various processing procedures, and necessary data, and executes various processes based on these programs. . The control unit 102 roughly includes a data acquisition unit 102a, an expression region estimation unit 102b, a gene structure prediction unit 102c, an expression value calculation unit 102d, and a prediction target region determination unit 102e.

  The data acquisition unit 102a uses a tiling array in which partial base sequences regularly extracted from the genomic base sequence are arranged as probes, and tiling array data relating to the expression intensity of each probe measured by the tiling array. And genome base sequence data related to the genome base sequence is acquired.

  The expression region estimation unit 102b estimates the gene expression region from the genome base sequence based on the tiling array data acquired by the data acquisition unit 102a. Specifically, as shown in FIG. 4, the expression region estimation unit 102b further includes an expression region candidate determination unit 102b1, an expression region candidate determination unit 102b2, and an expression region determination unit 102b3. The expression region candidate determination unit 102b1 determines expression region candidates that are gene expression region candidates based on the tiling array data. Specifically, as shown in FIG. 5, the expression region candidate determination unit 102b1 further includes a median value calculation unit 102b11, a region movement unit 102b12, and an expression region candidate selection unit 102b13. The median value calculation unit 102b11 calculates a median value of the expression intensities of probes included in a region having a predetermined length in the genome base sequence. The region moving unit 102b12 moves the region for which the median value is calculated by the median value calculating unit 102b11 along the genome base sequence. The expression region candidate selection unit 102b13 selects the maximum one from the plurality of median values accumulated by repeatedly executing the processing performed by the median value calculation unit 102b11 and the region movement unit 102b12, and selects the selected maximum median value. The region that was the calculation target is selected as an expression region candidate. Returning to FIG. 4 again, the expression region candidate determination unit 102b2 determines whether or not the expression region candidate determined by the expression region candidate determination unit 102b1 is statistically significant. When the expression region candidate determination unit 102b2 determines that the expression region candidate determination unit 102b3 is significant, the expression region determination unit 102b3 determines the expression region candidate as the expression region.

  Returning to FIG. 1 again, the gene structure prediction unit 102c uses the base corresponding to the expression region estimated in the expression region estimation step as a starting point, and adds tiling array data and genomic base sequence data for each base adjacent to the starting point. The gene structure is predicted by determining attributes (exons, introns, etc.) related to the structure of the gene. Specifically, as shown in FIG. 6, the gene structure prediction unit 102c further includes an expression region base probability determination unit 102c1, an adjacent base probability determination unit 102c2, a likelihood determination unit 102c3, and an attribute determination unit 102c4. . The expression region base probability determining unit 102c1 determines, using a Markov model, the probability related to the attribute of the expression region base for the expression region base that is the base corresponding to the expression region estimated by the expression region estimation unit 102b. The adjacent base probability determining unit 102c2 determines the state transition probability related to the attribute of the adjacent base using the Markov model for the adjacent base adjacent to the expression region base targeted by the expression region probability determining unit 102c1. The likelihood determining unit 102c3 repeatedly executes the adjacent base probability determining unit 102c2 up to the terminal base of the genome base sequence, and then determines the likelihood related to the attribute of each base using the maximum likelihood method. The attribute determination unit 102c4 combines the attribute of each base based on the probability determined by the expression region base probability determination unit 102c1, the state transition probability determined by the adjacent base probability determination unit 102c2, and the likelihood determined by the likelihood determination unit 102c3. Confirm.

  Returning to FIG. 1 again, the expression value calculation unit 102d calculates the expression value of the gene based on the expression intensity of the probe corresponding to the exon region included in the gene structure predicted by the gene structure prediction unit 102c.

  The prediction target region determination unit 102e determines, as a prediction target region, an unpredicted region that is a genomic base sequence excluding the region of the gene whose structure is predicted in the gene structure prediction step. Specifically, as illustrated in FIG. 7, the prediction target region determination unit 102e further includes an unpredicted region extraction unit 102e1, a base length determination unit 102e2, and a prediction target region determination unit 102e3. Here, the unpredicted region extraction unit 102e1 extracts the unpredicted region from the genome base sequence. The base length determination unit 102e2 calculates the base length of the unpredicted region extracted by the unpredicted region extraction unit 102e1, and determines whether the calculated base length is equal to or less than a predetermined base length. When the base length determination unit 102e2 determines that the prediction target region determination unit 102e3 is not less than or equal to the predetermined base length, the prediction target region determination unit 102e3 determines the unpredicted region extracted by the unpredicted region extraction unit 102e1 as the prediction target region.

  In the above configuration, main processing performed by the control unit 102 of the gene structure prediction apparatus 100 will be described with reference to FIGS. FIG. 8 is a flowchart showing an example of the main process.

  First, the data acquisition unit 102a acquires tiling array data and genome base sequence data, stores the acquired genome base sequence data in the base sequence data file 106a, and stores the acquired tiling array data in the tiling array data file 106b. (Step SA-1).

  Next, the expression region estimation unit 102b estimates the gene expression region from the genome base sequence based on the tiling array data acquired in step SA-1 (step SA-2: expression region estimation process). Specifically, the expression region estimation unit 102b determines expression region candidates that are gene expression region candidates based on the tiling array data acquired in step SA-1, and the determined expression region candidates are statistically determined. It is determined whether or not it is significant, and when it is determined that it is significant, an expression region candidate is determined as an expression region.

数式１において、ｖ_i（上付き線あり）は塩基番号ｉの位置において設定されたウインドウに含まれるプローブの発現強度ｖ_iの中央値である。また、ｋは、プライマリーポイント（ｐｒｉｍａｒｙ ｐｏｉｎｔ）と呼ばれるものであり、１以上Ｎ（ゲノム塩基配列の長さ）以下の整数値である。なお、ｍａｘｖ_i（上付き線あり）は推定の発現値（ｐｕｔａｔｉｖｅ ｅｘｐｒｅｓｓｉｏｎ ｖａｌｕｅ）と呼ばれるものである。 Here, the expression region estimation process performed by the expression region estimation unit 102b will be described with reference to FIG. FIG. 9 is a flowchart showing an example of the expression region estimation process.
First, the expression region candidate determination unit 102b1 determines expression region candidates that are gene expression region candidates based on the tiling array data (step SB-1). Specifically, expression region candidates are determined according to the following procedures (1) to (3).
(1) The median value calculation unit 102b11 calculates a median value of expression intensities of probes included in a region having a predetermined length (window: for example, 750 base length) in the genome base sequence.
(2) The region moving unit 102b12 moves the region for which the median value is calculated by the median value calculating unit 102b11 along the genome base sequence.
(3) In the expression region candidate selection unit 102b13, the median calculation unit 102b11 and the region movement unit 102b12 are repeatedly executed to select the maximum one from the accumulated median values, and the maximum median value selected The region that was the calculation target is selected as an expression region candidate.
That is, the primary point k shown in Formula 1 is obtained by the above procedure.

In Equation 1, v _i (with superscript line) is the median value of the expression intensity v _i of the probe included in the window set at the position of the base number i. Further, k is called a primary point and is an integer value of 1 or more and N (length of genome base sequence). Note that maxv _i (with a superscript line) is called an estimated expression value (putty expression value).

数式２において、Ｍ_kはプライマリーポイントｋに対応するウインドウに含まれるプローブの総数（例えば３０）である。また、ｍ_kは、プライマリーポイントｋに対応するウインドウに含まれるプローブのうち、発現強度が予め設定した閾値θより大きいプローブの数である。また、πθは、プローブの発現強度が閾値θより大きい確率であり、下記の数式３で定義される。

数式３において、δ_iは下記の数式４で定義される値である。また、λ_iは、１塩基あたりのプローブ数であり、下記の数式５で定義される。また、Ｉｆ［ｃｏｎｄｉｔｉｏｎ］は、括弧内が真であれば１を、括弧内が偽であれば０を返す関数である。また、Ｎはゲノム塩基配列の長さである。

数式４において、ｖ_iは塩基番号ｉの塩基に対応するプローブの発現強度である。θは予め定めた閾値である。なお、θの初期値には、全プローブの発現強度の中央値を設定してもよい。

数式５において、ｎ_iは塩基番号ｉの塩基に対応するプローブの数である。また、ｌはプローブの長さ（塩基数）である。
（３）数式２で算出した確率が期待値（Ｍ_k×πθ）より有意に大きい場合には、下記の数式６に示す有意水準で帰無仮説を棄却する。つまり、ステップＳＢ−１で決定した発現領域候補が統計的に有意であると判定する。これにより、プライマリーポイントｋに対応するウインドウに含まれるプローブの発現強度が偶然に観測されたものでないことが、統計的検定により示された。

数式６において、Ｍ_iは塩基番号ｉに対応するウインドウに含まれるプローブの総数である。また、ｍ_iは、塩基番号ｉに対応するウインドウに含まれるプローブのうち、発現強度が予め設定した閾値θより大きいプローブの数である。また、πθは、プローブの発現強度が閾値θより大きい確率であり、数式３で定義される。 Next, the expression region candidate determination unit 102b2 determines whether or not the expression region candidate determined by the expression region candidate determination unit 102b1 is statistically significant (step SB-2). Specifically, the determination is made according to the following procedures (1) to (3).
(1) A null hypothesis is set that the expression intensity of the probe included in the window corresponding to the primary point k determined in step SB-1 is observed by chance.
(2) The probability relating to the total number M _{k of} probes included in the window and the number m _k of probes whose expression intensity is greater than a preset threshold θ is calculated. Here, the probability follows a binomial distribution shown in Equation 2 below.

In Equation 2, M _k is the total number (for example, 30) of probes included in the window corresponding to the primary point k. M _k is the number of probes whose expression intensity is greater than a preset threshold θ among probes included in the window corresponding to the primary point k. Further, πθ is a probability that the expression intensity of the probe is larger than the threshold value θ, and is defined by Equation 3 below.

In Equation 3, δ _i is a value defined by Equation 4 below. Also, λ _i is the number of probes per base and is defined by Equation 5 below. If [condition] is a function that returns 1 if the parentheses are true and returns 0 if the parentheses are false. N is the length of the genome base sequence.

In Equation 4, v _i is the expression intensity of the probe corresponding to the base with the base number i. θ is a predetermined threshold value. The median value of the expression intensity of all probes may be set as the initial value of θ.

In Equation 5, n _i is the number of probes corresponding to the base with the base number i. L is the length (number of bases) of the probe.
(3) If the probability calculated by Equation 2 is significantly greater than the expected value (M _k × πθ), the null hypothesis is rejected at the significance level shown in Equation 6 below. That is, it is determined that the expression region candidate determined in step SB-1 is statistically significant. Thereby, it was shown by the statistical test that the expression intensity of the probe included in the window corresponding to the primary point k was not accidentally observed.

In Equation 6, M _i is the total number of probes included in the window corresponding to the base number i. M _i is the number of probes whose expression intensity is greater than a preset threshold θ among probes included in the window corresponding to the base number i. Further, πθ is a probability that the expression intensity of the probe is larger than the threshold θ, and is defined by Equation 3.

  Next, when the expression region determination unit 102b3 determines that the expression region is significant in Step SB-2, the expression region candidate determined in Step SB-1 is determined as an expression region (Step SB-3).

In the process of the gene structure prediction unit 102c, the gene structure (specifically, the exon region and the intron) is determined based on the assumption that the exon region corresponding to the primary point k is actually expressed. Region).
This is the end of the description of the expression region estimation process.

数式７において、Ｚ（ｓ，θ）は、下記の数式８で定義されるオッズである。また、Ｚ_T（ｓ）は、エクソン領域を尤もらしい長さに保つための式であり、下記の数式１５で定義される。また、Ｚ_U（ｓ）は、イントロン領域を尤もらしい長さに保つための式であり、下記の数式１８で定義される。

数式８において、Ｐ（ｘ，ｓ，δ｜θ）は、遺伝子である確率であり、下記の数式９で定義される、なお、Ｐ（ｘ，ｓ，δ｜θ）は、式の変形を行うことで最終的に下記の数式１０で表すことができる。また、Ｐ（ｘ，ｓ（上付き線あり），δ｜θ）は、遺伝子でない確率であり、式の変形を同様に行うことで最終的に下記の数式１１で表すことができる。

Returning to FIG. 8 again, the base structure corresponding to the expression region estimated in step SA-2 is used as a starting point in the gene structure prediction unit 102c, and the tiling array data and the genomic base sequence data are obtained for each base that is adjacent to the starting point. Based on the attribute (exon, intron, etc.) related to the structure of the gene, the structure of the gene is predicted (step SA-3: gene structure prediction process). Specifically, in the gene structure prediction unit 102c, (1) using the Markov model for the expression region base that is the base corresponding to the expression region estimated in step SA-2, the probability relating to the attribute of the expression region base is used. (2) For the adjacent bases adjacent to the expression region base, the state transition probability relating to the attribute of the adjacent base is determined using a Markov model, and the processing of (3) (2) is performed as a genomic base sequence After repeatedly executing up to the base of the base, the likelihood regarding the attribute of each base is determined using the maximum likelihood method, and (4) the combination of the attributes of each base based on the determined probability, state transition probability and likelihood Confirm. More specifically, the gene structure prediction unit 102c calculates the score function shown in Equation 7 below by dynamic programming (or may be Newton's method or Powell's method), thereby maximizing the value of the score function. And the threshold θ are determined to predict the gene structure. Here, the attribute s corresponds to “state number” in FIG.

In Equation 7, Z (s, θ) is an odds defined by Equation 8 below. Z _T (s) is an expression for keeping the exon region at a reasonable length, and is defined by Expression 15 below. Z _U (s) is an expression for maintaining the intron region to have a reasonable length, and is defined by Expression 18 below.

In Equation 8, P (x, s, δ | θ) is a probability of being a gene, and is defined by Equation 9 below, where P (x, s, δ | θ) By doing so, it can be finally expressed by the following formula 10. Further, P (x, s (with superscript line), δ | θ) is a probability that it is not a gene, and can be finally expressed by the following formula 11 by similarly modifying the formula.

Here, the gene structure prediction process performed by the gene structure prediction unit 102c will be described with reference to FIG. FIG. 10 is a flowchart showing an example of gene structure prediction processing.
First, the expression region base probability determination unit 102c1 determines, using a Markov model, the probability related to the attribute of the expression region base for the expression region base that is the base corresponding to the expression region estimated in step SA-2. Step SC-1). Specifically, for the primary point k, the attribute s _k is set for the probabilities P (x _k , s _k ) and P (x _k , s _k (with superscript lines)) included in the equations 10 and 11. While changing, the maximum value of “P (x _k , s _k ) ÷ P (x _k , s _k (with superscript line))” is determined.

Next, the adjacent base probability determination unit 102c2 determines the state transition probability related to the attribute of the adjacent base using the Markov model for the adjacent base adjacent to the expression region base targeted in step SC-1 (step SC). -2: See FIG. Specifically, the state transition probabilities P (x _{k + 1} , s _{k + 1} | x _k , s _k ) and P (x _{k + 1} , s _{k + 1)} (superscript lines included in Equations 10 and 11) Yes) | x _k , s _k (with superscript line)), while changing the attribute s _{k + 1} , “P (x _{k + 1} , s _{k + 1} | x _k , s _k ) ÷ P (x _{k +1} , s _{k + 1} (with superscript line) | x _k , s _k (with superscript line)) ”is determined. Similarly, the state transition probabilities P (x _k−1 , s _k−1 | x _k , s _k ) and P (x _k−1 , s _k−1 (superscript lines) included in the equations 10 and 11 are used. Yes) | x _k , s _k (with superscript line)), while changing the attribute s _k-1 , “P (x _k−1 , s _k-1 | x _k , s _k ) ÷ P (x _{k −1} , s _k-1 (with superscript line) | x _k , s _k (with superscript line)) ”is determined. Here, in FIG. 12, 3 ′ end and 5 ′ end are the end points of the state transition.

Then, the process of step SC-2 is repeated until the terminal base of the genome base sequence. As a result, the product of the state transition probabilities P (x _i , s _i | x _{i + 1} , s _{i + 1} ) included in the equations 10 and 11 and the state transition probabilities P (x _i , s _i (with superscript lines) ) | X _{i + 1} , s _{i + 1} (with superscript line)), the product of “P (x _i , s _i | x _{i + 1} , s _{i + 1} ) ÷ P (x _i , s _{i ”} (with superscript line) | x _{i + 1} , s _{i + 1} (with superscript line))” is determined. Similarly, the state transition probability _{P (x i, s i |} x i-1, s i-1) included in Equation 10 and Equation 11 the product and the state transition probability P (x _i, s _i (superscript With the line) | x _i-1 , s _i-1 (with superscript line)), the product of “P (x _i , s _i | x _i-1 , s _i-1 ) ÷ P (x _i , S _i (with superscript line) | x _i−1 , s _i-1 (with superscript line) ”is determined. In FIG. 12, the extension of the region from the start point to the end point is, for example, at each position over a continuous base sequence in which the end point state (the state of base number 0 and base number 25) is length L (for example, 800 base length). Stop if you keep the maximum score at.

ここで、各遺伝子は、異なる発現レベルをもつので、閾値θおよびＰ（δ_i｜ｓ_i，θ）は、各転写単位（ＴＵ：Ｔｒａｎｓｌａｔｉｏｎａｌ Ｕｎｉｔ）内で局所的に決めるべきである。従って、数式１２において、ａ、ｂ、ｃ、およびｄは、それぞれ下記の数式１３で定義される。

数式１３において、ａ、ｂ、ｃ、およびｄは、５'末端および３'末端間の領域内において決められる。

Next, the likelihood determining unit 102c3 determines the likelihood related to the attribute of each base using the maximum likelihood method (step SC-3). Specifically, regarding the product of likelihood P (δ _i | s _i , θ) and the product of likelihood P (δ _i | s _i (with superscript), θ) included in Equation 10 and Equation 11, The value of “product of P (δ _i | s _i , θ) ÷ product of P (δ _i | s _i (with superscript), θ)” is determined. For example, in FIG. 12, the value of P (δ _i | s _i , θ) is updated for the region from base number 0 to base number 25. Here, P (δ _i | s _i , θ) and P (δ _i | s _i (with superscript line), θ) are defined by the following equations 12 and 14, respectively.

Here, since each gene has a different expression level, the threshold values θ and P (δ _i | s _i , θ) should be determined locally within each transcription unit (TU). Therefore, in Expression 12, a, b, c, and d are respectively defined by Expression 13 below.

In Equation 13, a, b, c, and d are determined in the region between the 5 ′ end and the 3 ′ end.

  With the above processing, the maximum value of Z (s, θ) included in Equation 7 is determined at the preset threshold θ.

数式１５において、Ｎ_Tはエクソン数である。また、ｌ_T,jはエクソン長である。また、Ｐ_T（ｌ）はエクソン長ｌに関する確率分布であり、下記の数式１６で定義される。また、Ｅ_TはＰ_T（ｌ）の期待値であり、下記の数式１７で定義される。

数式１６において、μ_Tはエクソン長の平均値である。σ_Tはエクソン長の標準偏差である。なお、μ_Tおよびσ_Tはトレーニングセットで観測された実測値である。

数式１８において、Ｎ_Uはイントロン数である。また、ｌ_U,jはイントロン長である。また、Ｐ_U（ｌ）はイントロン長ｌに関する確率分布であり、下記の数式１９で定義される。また、Ｅ_UはＰ_U（ｌ）の期待値であり、下記の数式２０で定義される。

数式１９において、μ_Uはイントロン長の平均値である。σ_Uはイントロン長の標準偏差である。なお、μ_Uおよびσ_Uはトレーニングセットで観測された実測値である。

Then, the value of the expression Z _T (s) for keeping the exon region shown in Equation 15 below to have a reasonable length and the equation Z _U (s) for keeping the intron region shown in Equation 18 below in a reasonable length. ) Is finally determined, the value (maximum value) of the score function of Equation 7 at the preset threshold θ is determined, and the maximum value of the determined score function, the threshold θ and the maximum value are determined. The obtained attributes s (s ₁ ,... S _k ,... S _N ) are associated with each other and stored in a predetermined area of the gene structure data file 106d.

In Equation 15, N _T is the number of exons. L _{T, j} is the length of the exon. P _T (l) is a probability distribution with respect to the exon length l and is defined by Equation 16 below. E _T is an expected value of P _T (l) and is defined by the following Equation 17.

In Equation 16, μ _T is an average value of exon lengths. σ _T is the standard deviation of the exon length. Μ _T and σ _T are actually measured values observed in the training set.

In Equation 18, N _U is the number of introns. L _{U, j} is the intron length. P _U (l) is a probability distribution with respect to the intron length l and is defined by Equation 19 below. E _U is an expected value of P _U (l), and is defined by the following Equation 20.

In Equation 19, μ _U is the average value of the intron length. σ _U is the standard deviation of the intron length. Μ _U and σ _U are actually measured values observed in the training set.

  Next, the threshold θ is reset to another value, and the above processing is performed in the same manner.

Next, the attribute determination unit 102c14 determines the attribute combination of each base based on the probability, state transition probability, and likelihood determined in step SC-3 (step SC-4). Specifically, the maximum value of the score function at each threshold value θ accumulated in the gene structure data file 106d is determined, and the attribute s corresponding to the determined maximum value is determined as the optimum value. As a result, the combination of the threshold value θ and the attribute s that maximizes the score function value was determined, and the gene structure was predicted.
This is the end of the description of the gene structure prediction process.

  Returning to FIG. 8 again, the prediction target region determination unit 102e determines an unpredicted region, which is a genomic base sequence excluding the region of the gene whose structure is predicted in step SA-3, as a prediction target region (step SA-4). : Prediction target area determination processing). Specifically, the prediction target region determination unit 102e extracts an unpredicted region from the genome base sequence, calculates the base length of the extracted unpredicted region, and whether the calculated base length is equal to or less than a predetermined base length. If it is determined that the length is not less than the predetermined base length, the extracted unpredicted region is determined as the prediction target region.

Here, the prediction target region determination process performed by the prediction target region determination unit 102e will be described with reference to FIG. FIG. 11 is a flowchart illustrating an example of the prediction target area determination process.
First, the unpredicted region extraction unit 102e1 extracts the unpredicted region from the genome base sequence (step SD-1).
Next, the base length determination unit 102e2 calculates the base length of the unpredicted region extracted in step SD-1, and determines whether the calculated base length is a predetermined base length (for example, 1000 base lengths) or less. (Step SD-2).
Next, when the prediction target region determination unit 102e3 determines that the length is not less than the predetermined base length in Step SD-2 (Step SD-3: Yes), the unpredicted region extracted in Step SD-1 is determined as the prediction target region. (Step SD-4). When it is determined in step SD-2 that the length is equal to or shorter than the predetermined base length (step SD-3: No), it is determined that there cannot be a gene structure that falls within the region of the predetermined base length, and extracted. The unpredicted area is treated as an area whose structure is predicted, and the process returns to step SD-1.
This is the end of the description of the prediction target region determination process.

  Returning to FIG. 8 again, when the prediction target area is determined in step SA-4 (step SA-5: No), the process returns to step SA-2 again, and the prediction target area is not determined in step SA-4. If this is the case (step SA-5: Yes), the main process is terminated.

  Here, the expression value calculation unit 102d calculates the expression value of the gene based on the expression intensity of the probe corresponding to the exon region included in the structure of the gene predicted in step SA-3, and relates to the calculated expression value. The data may be stored in a predetermined storage area of the expression value data file 106e.

  This completes the description of the main process.

  As described above, the gene structure prediction apparatus 100 predicts a gene structure from a genome base sequence based on tiling array data. Specifically, the gene structure prediction apparatus 100 acquires tiling array data and genomic base sequence data, estimates a gene expression region from the genomic base sequence based on the acquired tiling array data, and estimates The structure of the gene is predicted by determining the attributes related to the structure of the gene based on the tiling array data and the genomic base sequence data for each base that is adjacent to the base corresponding to the expressed region. Therefore, the region and structure of an unknown gene can be accurately predicted from the genome base sequence based on the tiling array data. The present invention can also be applied to non-coding RNA. Here, as shown in FIG. 14, the exons and introns in the predicted gene structure were highly correlated with the true exons and introns based on the cDNA information. And the correlation coefficient was higher than the correlation coefficient in the conventional method which predicts a gene structure by determining whether it is an exon or an intron depending on whether a threshold value is exceeded. Furthermore, as shown in FIG. 15, the predicted gene structure (the gene structure indicated by MA-3-2 in FIG. 15 (orange represents an exon region and light orange represents an intron region)). It substantially corresponded to the true gene structure (the gene structure shown in MA-3-1 in FIG. 15 (blue represents an exon region and light blue represents an intron region)).

  Here, FIG. 15, which is a display screen displaying the predicted gene structure and the true gene structure so as to be comparable, will be briefly described. In FIG. 15, MA-1 is a display area that displays an illustration of the target genome, a predicted gene region (8884876 bp to 8865453 bp: 16557 bp in MA-1), and the like. MA-2 is an operation button for enlarging / reducing a gene structure displayed in a display area MA-3 described later. MA-3 is a display area for displaying the gene structure along the genome base sequence. In particular, MA-3-1 is a region that displays the true gene structure to be compared, and MA-3-2 is a region that displays the predicted gene structure along the base sequence. In MA-3-1, the exon region included in the true gene structure is represented in blue, and the intron region is represented in light blue. In MA-3-2, exon regions included in the displayed gene structure are shown in orange, and intron regions are shown in light orange. Further, in MA-3-2, along the base sequence, the expression intensity of the probe corresponding to each base is displayed in green when it is equal to or lower than the median, and is displayed in red when it exceeds the median. MA-3-3 is an operation button used for switching data to be displayed.

  In the gene structure prediction apparatus 100, the expression region estimation unit 102b determines expression region candidates that are gene expression region candidates based on the tiling array data, and the determined expression region candidates are statistically significant. If the expression region candidate is determined to be significant, the expression region candidate is determined as the expression region. Therefore, the expression region can be estimated easily and accurately using an existing statistical method.

  Further, in the gene structure prediction apparatus 100, the expression region candidate determination unit 102b1 calculates the median of the expression intensity of the probes included in the region for a region of a predetermined length in the genome base sequence, By moving the region to be calculated along the genome base sequence and repeatedly executing these processes, the largest one is selected from the accumulated medians, and the selected median of the maximum median is selected. Since the selected region is selected as an expression region candidate, an appropriate expression region candidate can be selected in consideration of variations in expression intensity of each probe measured by the tiling array.

  Further, in the gene structure prediction apparatus 100, the gene structure prediction unit 102c determines, using a Markov model, a probability relating to an attribute of the expression region base for an expression region base that is a base corresponding to the estimated expression region, For neighboring bases adjacent to the target expression region base, the state transition probability related to the attribute of the neighboring base is determined using a Markov model, and the state transition probability is repeatedly determined up to the end base of the genome base sequence. After that, the likelihood regarding the attribute of each base is determined using the maximum likelihood method, and the combination of the attributes of each base is determined based on the determined probability, state transition probability, and likelihood. Further, it can be predicted with high accuracy.

  Moreover, since the gene structure prediction apparatus 100 calculates the expression value of the gene based on the expression intensity of the probe corresponding to the exon region included in the predicted gene structure, the expression level of the unknown gene can be obtained. it can.

  Further, the gene structure prediction apparatus 100 determines an unpredicted region, which is a genomic base sequence excluding a region of a gene whose structure is predicted, as a prediction target region, and the above-described expression region estimation for the determined prediction target region Since the process and the gene structure prediction process are executed again, the regions and structures of a plurality of genes can be efficiently predicted from the genome base sequence.

  In the gene structure prediction apparatus 100, the prediction target region determination unit 102e extracts an unpredicted region from the genome base sequence, calculates the base length of the extracted unpredicted region, and the calculated base length is equal to or less than a predetermined base length. When it is determined whether or not it is not less than a predetermined base length, the extracted unpredicted region is determined as a prediction target region, so that a short region where the gene structure cannot enter can be excluded, As a result, a prediction target area suitable for the structure prediction can be determined.

Finally, a procedure for transforming the above-described Expression 9 into Expression 10 will be described. In Equation 9, P _k is defined by Equation 21 below, P _k → _N is defined by Equation 22 below, and P ₁ → _k is defined by Equation 23 below.

First, Equation 21, Equation 22, and Equation 23 can be transformed into Equation 24, Equation 25, and Equation 26, respectively, according to the law of posterior probability.

Next, since it is natural to assume that δ _i depends on s _i and θ and is not affected by x _i , in equations 24, 25 and 26, P (δ _k | x _k , s _k , θ) can be transformed into Equation 27 below.

Further, using the Markov model basic hypothesis that the probabilities of x _i and s _i in Equation 24 depend on adjacent x _i-1 and s _i-1 , P (x _k , s _k | θ in Equation 24 ) Can be transformed into Equation 28 below.

Similarly, in Equations 25 and 26, P (x _i , s _i | x _i−1 , s _i−1 , δ _i−1 , θ) and P (x _i , s _i | x _{i + 1} , s _{i + 1} , δ _{i + 1} , θ) can be transformed into the following equations 29 and 30, respectively.

Then, P _k defined by Equation 21 can be finally transformed into Equation 31 below.

Similarly, P _k → _N defined by Expression 22 can be finally transformed into the following Expression 32, and P ₁ → _k defined by Expression 23 can be finally transformed into the following Expression 33.

  Equation 10 can be derived by modifying Equation 9 in the above procedure.

  As described above, the gene structure prediction method and gene structure prediction program according to the present invention can accurately predict the region and structure of an unknown gene from a genomic base sequence based on tiling array data. It is extremely useful in fields such as medicine.

1 is a block diagram showing a configuration of a gene structure prediction apparatus 100. FIG. It is a figure which shows an example of the information stored in the base sequence data file 106a. It is a figure which shows an example of the information stored in the tiling array data file 106b. It is a block diagram which shows the structure of the expression area estimation part 102b. It is a block diagram which shows the structure of the expression area candidate determination part 102b1. It is a block diagram which shows the structure of the gene structure estimation part 102c. It is a block diagram which shows the structure of the prediction object area | region determination part 102e. It is a flowchart which shows an example of a main process. It is a flowchart which shows an example of an expression area estimation process. It is a flowchart which shows an example of a gene structure prediction process. It is a flowchart which shows an example of a prediction object area | region determination process. It is a figure which shows an example of the Markov model regarding arrangement | sequence data. It is a figure which shows the definition of the attribute s. It is a figure which shows an example of the evaluation result of prediction accuracy. It is a figure which shows an example of the estimated gene structure.

Explanation of symbols

100 Gene structure prediction device
102 Control unit
102a Data acquisition unit
102b Expression region estimation unit
102b1 expression region candidate determination unit
102b11 median value calculator
102b12 area moving part
102b13 Expression region candidate selection section
102b2 Expression region candidate determination unit
102b3 expression region determining unit
102c Gene structure prediction unit
102c1 expression region base probability determination unit
102c2 Adjacent base probability determination unit
102c3 likelihood determination unit
102c4 attribute determination part
102d Expression value calculation unit
102e Prediction area determination unit
102e1 unpredicted region extraction unit
102e2 base length determination unit
102e3 prediction target region determination unit
104 Communication interface
106 Storage unit
106a Genome sequence data file
106b Tiling array data file
106c Expression region data file
106d gene structure data file
106e Expression value data file
108 Input / output interface
112 Input device
114 Output device 200 External system 300 Network

Claims (10)

Using a tiling array in which partial base sequences regularly extracted from genomic base sequences are arranged as probes, tiling array data relating to the expression intensity of each probe measured with the tiling array and genome relating to the genomic base sequence A data acquisition step for acquiring base sequence data;
Based on the tiling array data acquired in the data acquisition step, an expression region estimation step for estimating a gene expression region from the genomic base sequence;
Starting from the base corresponding to the expression region estimated in the expression region estimation step, the attribute relating to the structure of the gene is determined for each base that is adjacent to the origin on the basis of the tiling array data and the genomic base sequence data. A gene structure prediction step for predicting the structure of the gene,
Only including,
The expression region estimation step includes:
An expression region candidate determination step for determining an expression region candidate which is a candidate for a gene expression region based on the tiling array data;
Expression region candidate determination step for determining whether the expression region candidate determined in the expression region candidate determination step is statistically significant,
When it is determined to be significant in the expression region candidate determination step, an expression region determination step for determining an expression region candidate as the expression region;
Further including
The expression region candidate determination step includes
For a region of a predetermined length in the genomic base sequence, a median calculation step for calculating the median of the expression intensity of the probes contained in the region;
A region moving step of moving the region that is the median calculation target in the median calculation step along the genome base sequence;
The median calculation step and the region movement step are repeatedly executed to select the largest one among the plurality of median values accumulated, and the region that was the target of calculation of the selected maximum median value is selected as the expression region candidate An expression region candidate selection step selected as:
Gene structure prediction method comprising the further contains Mukoto. The gene structure prediction step includes:
For the expression region base that is the base corresponding to the expression region estimated in the expression region estimation step, the expression region base probability determination step for determining the probability regarding the attribute of the expression region base using a Markov model,
For adjacent bases adjacent to the expression region base targeted in the expression region probability determination step, adjacent base probability determination step for determining a state transition probability related to the attribute of the adjacent base using a Markov model,
A likelihood determination step of determining the likelihood related to the attribute of each base using a maximum likelihood method after repeatedly executing the adjacent base probability determination step up to the terminal base of the genomic base sequence,
An attribute determination step for determining a combination of attributes of each base based on the probability, the state transition probability and the likelihood;
The gene structure prediction method according to claim 1 , further comprising: Based on the expression intensity of the probe corresponding to the exon region included in the gene structure predicted in the gene structure prediction step, an expression value calculation step for calculating the expression value of the gene,
Gene structure prediction method according to claim 1 or 2, further comprising a. A prediction target region determination step for determining an unpredicted region that is a genomic base sequence excluding a region of the gene whose structure is predicted in the gene structure prediction step, as a prediction target region;
Further including
Re-execution of the expression region estimation step and the gene structure prediction step on the prediction target region determined in the prediction target region determination step;
The gene structure prediction method according to any one of claims 1 to 3 , wherein: The prediction target region determination step includes
An unpredicted region extraction step of extracting the unpredicted region from the genome base sequence;
Calculating the base length of the unpredicted region extracted in the unpredicted region extraction step, and determining whether the calculated base length is equal to or less than a predetermined base length; and
When it is determined that the base length determination step is not less than or equal to a predetermined base length, a prediction target region determination step for determining the unpredicted region extracted in the unpredicted region extraction step as the prediction target region;
The gene structure prediction method according to claim 4 , further comprising: Using a tiling array in which partial base sequences regularly extracted from genomic base sequences are arranged as probes, tiling array data relating to the expression intensity of each probe measured with the tiling array and genome relating to the genomic base sequence A data acquisition step for acquiring base sequence data;
Based on the tiling array data acquired in the data acquisition step, an expression region estimation step for estimating a gene expression region from the genomic base sequence;
Starting from the base corresponding to the expression region estimated in the expression region estimation step, the attribute related to the structure of the gene is determined for each base that is adjacent to the origin based on the tiling array data and the genome base sequence data. A gene structure prediction step for predicting the structure of the gene,
Only including,
The expression region estimation step includes:
An expression region candidate determination step for determining an expression region candidate which is a candidate for a gene expression region based on the tiling array data;
Expression region candidate determination step for determining whether the expression region candidate determined in the expression region candidate determination step is statistically significant,
When it is determined to be significant in the expression region candidate determination step, an expression region determination step for determining an expression region candidate as the expression region;
Further including
The expression region candidate determination step includes
For a region of a predetermined length in the genomic base sequence, a median calculation step for calculating the median of the expression intensity of the probes contained in the region;
A region moving step of moving the region that is the median calculation target in the median calculation step along the genome base sequence;
The median calculation step and the region movement step are repeatedly executed to select the largest one among the plurality of median values accumulated, and the region that was the target of calculation of the selected maximum median value is selected as the expression region candidate An expression region candidate selection step selected as:
Further gene structure prediction program characterized by executing the including gene structure prediction methodologies to computers. The gene structure prediction step includes:
For the expression region base that is the base corresponding to the expression region estimated in the expression region estimation step, the expression region base probability determination step for determining the probability regarding the attribute of the expression region base using a Markov model,
For adjacent bases adjacent to the expression region base targeted in the expression region probability determination step, adjacent base probability determination step for determining a state transition probability related to the attribute of the adjacent base using a Markov model,
A likelihood determination step of determining the likelihood related to the attribute of each base using a maximum likelihood method after repeatedly executing the adjacent base probability determination step up to the terminal base of the genomic base sequence,
An attribute determination step for determining a combination of attributes of each base based on the probability, the state transition probability and the likelihood;
The gene structure prediction program according to claim 6 , further comprising: Based on the expression intensity of the probe corresponding to the exon region included in the gene structure predicted in the gene structure prediction step, an expression value calculation step for calculating the expression value of the gene,
The gene structure prediction program according to claim 6 or 7 , further comprising: A prediction target region determination step for determining an unpredicted region that is a genomic base sequence excluding a region of the gene whose structure is predicted in the gene structure prediction step, as a prediction target region;
Further including
Re-execution of the expression region estimation step and the gene structure prediction step on the prediction target region determined in the prediction target region determination step;
The gene structure prediction program according to any one of claims 6 to 8 , wherein: The prediction target region determination step includes
An unpredicted region extraction step of extracting the unpredicted region from the genome base sequence;
Calculating the base length of the unpredicted region extracted in the unpredicted region extraction step, and determining whether the calculated base length is equal to or less than a predetermined base length; and
When it is determined that the base length determination step is not less than or equal to a predetermined base length, a prediction target region determination step for determining the unpredicted region extracted in the unpredicted region extraction step as the prediction target region;
The gene structure prediction program according to claim 9 , further comprising:

Images