JP2007312653A - Analyzing method for character extraction and comparison classification of sequential gene expression data and analyzing apparatus based on the analyzing method - Google Patents

Abstract

<P>PROBLEM TO BE SOLVED: To provide an analyzing method for the extraction of characters of sequential gene expression data of a recombinant gene in individual plant body in the case of using a Luc-Tag line transformed plant and for the comparison classification based on the characters of the sequential gene expression data between the individuals, and provide an analyzing apparatus based on the analyzing method. <P>SOLUTION: Cluster analysis is carried out on multiple specimens in a line based on the wave form expressing the change of luciferase activity of a Luc-Tag line transformed plant with time, an individual group having high similarity is selected, a model waveform reflecting the similarity is formed, the model waveform is decomposed to a plurality of single-peaked waveforms and the characters of the change of luciferase activity with time in the individual group having high line similarity are grasped. <P>COPYRIGHT: (C)2008,JPO&INPIT

Description

  The present invention relates to feature extraction of time-series gene expression level data in individual organisms, analysis method for the purpose of comparative classification based on characteristics of time-series gene expression level data between individuals, and analysis apparatus based on the analysis method About. In particular, for transgenic transgenic plants, feature extraction of time-series gene expression data of the recombinant gene in individual plants and comparison based on characteristics of time-series gene expression data between individuals The present invention relates to an analysis method for classification and an analysis apparatus based on the analysis method.

  The temporal change in the expression level of a gene in an individual organism reflects the mechanism that controls the expression of the gene or the role of the gene product in the life activity of the organism. For example, in a series of steps of cell division called the so-called cell cycle, genes that code for specific types of proteins involved in the progress of that step are expressed in accordance with a pre-programmed time table. .

  In addition, higher organisms composed of multiple cells have various organs and produce specific gene products in the cells accompanying the process of forming individual organs, that is, the process of differentiation. is doing. For this reason, there have been reports of cases where the expression level of each specific gene in a differentiated cell changes with time depending on the stage reached by the cell. In higher organisms composed of multiple cells, the expression process of a gene encoding a protein present on genomic DNA includes, for example, transcription, maturation (splicing) into mRNA, translation based on mRNA, Through the steps of folding the translated peptide chain and post-translation processing, the target mature protein is produced. Each of these series of steps is a reaction involving an enzyme protein, and is influenced by the enzyme reaction rate. In general, the enzyme activity of an enzyme protein has temperature dependence, and the enzyme reaction rate varies depending on the temperature in the cell.

  For example, in a plant cell that is cultured at a predetermined temperature, the temperature in the cell is kept constant, and thus the enzyme activity of the enzyme protein itself is generally maintained at a constant level. Yes. On the other hand, in a plant body in an environment where the ambient temperature varies with time, the temperature in each cell of the plant body is affected by the ambient temperature and similarly varies with time. Correspondingly, the enzyme activity also shows temporal variations. Moreover, in the plant body which has photosynthesis ability, ATP and NADPH which are produced using light energy in its own photosynthesis process are used for a considerable part of the energy source consumed in order to maintain the life activity of the plant body. We are using. For example, in the synthesis of various amino acids and ribonucleoside 5 ′ triphosphates, which are substrate substances used in the expression process of a gene encoding a protein, the photosynthetic mechanism in the cell of a plant having photosynthesis ability ATP and NADPH produced in Japan are used in a substantial part. Therefore, when the amount of light irradiation with respect to the plant body varies with time, the amounts of various amino acids and ribonucleoside 5 'triphosphate synthesized in the cell also vary with time. For example, the life activity of a plant having photosynthetic ability, where the amount of light irradiation and the ambient temperature are periodically changing in synchronization with time, also show corresponding periodic fluctuations. Is predicted. In a cultivated environment that undergoes daytime and nighttime temperature fluctuations under sunlight irradiation, for example, in a plant that is cultivated, it seems to synchronize with an extrinsic cycle, including the effects of periodic variations due to the above factors, for example. There are also cases where a temporal change in the expression level of the gene is induced.

  In addition, in cells constituting various tissues and organs of a plant body, transcription from a specific gene itself is often induced by some transcription factor. For example, when an external stimulus induces an increase or activation of a transcription factor, and as a result, transcription from the gene is induced, while when an external stimulus is removed, an increase in the transcription factor thereafter, Activation does not proceed, and the amount of transcription often decreases with time. Also in this case, a temporal change in the expression level of the target gene occurs.

  Many genes encoding different proteins exist on the genomic DNA of higher organisms composed of multiple cells. At that time, the expression level of each gene is considered to exhibit a unique temporal change depending on the biological function of the protein expressed from those genes and the mechanism of the function expression. In other words, when studying the biological functions and roles of proteins, which are products derived from each gene, information on temporal changes in the expression level of the gene is also important information.

  On the other hand, in higher organisms composed of multiple cells, the protein produced by the expression of the coding gene in the cell once differentiated, with the passage of time, depending on its biological function, It loses its physiological activity (inactivation) and is further degraded by endogenous proteolytic enzymes. The concentration of the protein present in the cell along with the protein metabolism mechanism, which involves protein synthesis and degradation processes in the cell, shows temporal fluctuations. For example, when the expression of a specific gene is promoted by some factor in the cell, and the factor that promotes the expression of the gene is removed after the concentration of the protein encoded by the gene is increased, Due to protein degradation, the protein concentration decreases. As a result, the intracellular concentration of the protein shows a temporal change in the peak shape showing the maximum at a certain time.

After transcription and maturation (splicing) into mRNA, mRNA used for translation of the peptide chain of the protein is also used for translation and then degraded by endogenous nuclease (RNase). The concentration of the mRNA present in the cell in accordance with the metabolic mechanism of the mRNA, in which the production and degradation processes of the mRNA are involved, shows a temporal variation. For example, when the expression of a specific gene is promoted by some factor in the cell and the concentration of the mRNA is increased, and then the factor that promotes the expression of the gene is removed, the enzyme is subsequently degraded enzymatically. Therefore, the concentration of the mRNA decreases. As a result, the intracellular concentration of the mRNA shows a temporal change in the peak shape showing the maximum at a certain time. In general, the enzymatic degradation of mRNA proceeds much more rapidly than the enzymatic degradation of a protein composed of a peptide chain translated therefrom, and therefore the time at which the intracellular concentration of the mRNA reaches its maximum corresponds. It is slightly earlier than the time when the intracellular concentration of the protein shows the maximum. In addition, after the time showing the maximum, the intracellular concentration of the mRNA decreases rapidly.

When the intracellular concentration of a specific protein shows a temporal change in the peak shape as described above, higher organisms composed of multiple cells, for example, various tissues and organs possessed by a plant have a plurality of the same type of cells. However, when the cells of the same type are compared with each other, the change in the intracellular concentration of the protein is not completely synchronized between cells. However, considering the whole population of allogeneic cells contained in these specific tissues and organs, that is, there is still an “average value” that averages the intracellular concentration of the protein in the same type of individual cells. It has a tendency to show a temporal change in the peak shape indicating the maximum at time. That is, the entire population of the same type of cells contained in a specific tissue or organ has a tendency that the time at which the intracellular concentration of the protein shows a peak in each cell is concentrated within a specific time width. In this case, the “average value” of the entire population of allogeneic cells shows a maximum during this particular time span.

  For example, in a plant body, it is technically difficult to trace changes in the intracellular concentration of a target protein for individual cells constituting various tissues and organs of the plant individual. For this reason, the information derived by tracking the temporal change of the “average value” that averages the intracellular concentration of the protein in the same type of cell population contained in these specific tissues and organs is the product derived from each gene. It is used to study biological functions and roles of proteins.

  In studying the biological function and role of a protein, which is a product derived from each gene, it is most desirable to track temporal changes in the expression level of the gene in individual cells. However, it is technically difficult to trace temporal changes in the expression level of a specific gene in individual cells constituting each individual in a higher organism composed of multiple cells. Therefore, among a large number of cells constituting each individual, for example, an entire population of cells of the same kind that constitute various tissues and organs by differentiation, the expression of the gene in each cell included in the population A method of tracking the temporal change of the “average value” obtained by averaging the quantities is used.

  For example, for a whole group of cells of the same type that constitute various tissues and organs of the plant individual while continuing the growth of the plant body, a specific gene in each cell included in the group over time In order to track an “average value” obtained by averaging the expression level of a gene, conventionally, a part of cells is collected from the population, and the concentration of mRNA contained in the collected cell sample or a corresponding value is collected. Means for analyzing protein concentration have been utilized. This method of sampling a part of the cells is applicable only when the premise that the whole population of the same kind of cells constituting the population is in a somewhat homogeneous state is satisfied. Moreover, it is an effective means when the sampling operation does not act on the population after sampling a part of cells as an external stress.

  In fact, for a growing plant body, from a whole population of cells of the same kind contained in the plant individual, at the same time, a plurality of collected cell samples are collected, and the concentration of mRNA in each sample, or Comparing the results of analysis of the corresponding protein concentration, it was found that the “average value” of the intracellular concentration of mRNA and the “average value” of the intracellular concentration of the corresponding protein showed considerable dispersion in each sample. It turns out that there are not a few. That is, for the growing plant body, the “average value” of the intracellular concentration of mRNA and the “average value” of the intracellular concentration of the corresponding protein in the entire population of cells of the same kind contained in the plant individual. When tracking temporal changes in a time series over a considerably long time, it is determined that the technique of sampling some cells is not effective.

  Instead of the “destructive evaluation method” in which some cells are sampled, the present inventors and their collaborators are completely synchronized with the expression of the gene of interest as a non-destructive evaluation method, A luciferase gene that encodes a photoprotein (luciferase) is introduced so that the expression level is quantitatively proportional to the expression level, and the recombinant expression type photoprotein (luciferase) is produced in the cell, and the recombinant expression It was demonstrated that a technique for measuring chemiluminescence derived from a type of photoprotein can be employed.

  Specifically, for a specific gene present on the genomic DNA of a plant, a sequence encoding a protein amino acid sequence encoded by the gene (IRES) that continues transcription of pre-mRNA Immediately after, the region (ORF portion) encoding the photoprotein (luciferase) of the luciferase gene is inserted and ligated. This chimeric gene has a promoter sequence that causes transcription of a specific gene in the upstream portion, and transcription is similarly initiated by an endogenous factor that induces expression of the specific gene. In that case, the obtained mRNA becomes translatable into a peptide chain of a photoprotein (luciferase) encoded by the inserted luciferase gene. That is, translation is performed based on the mRNA, and when folding proceeds into a protein, an active photoprotein (luciferase) is produced in the cell as a reporter protein. The expression level corresponds to the expression level of the target protein encoded by the target specific gene in the wild-type plant.

The concentration of the reporter protein (luciferase) present in the cells can be quantitatively determined by evaluating the enzyme activity of the photoprotein (luciferase). Specifically, when the substrate substance of the photoprotein (luciferase) is present in the cell (l) at a predetermined concentration C _sub. , The intracellular concentration C _fusion (l) of the fusion protein indicates that C terminal The chemiluminescence per unit time due to the enzyme activity of a partial photoprotein (luciferase): dP _chem. (L) / dt is calculated using the apparent rate constant k _react.
dP _chem. (l) / dt = k _react. C _sub. C _fusion (l) (1)
Can be displayed. Therefore, chemiluminescence per unit time: dPchem. By measuring / dt, the intracellular concentration C _fusion (l) of the fusion protein can be estimated. When targeting the entire cell population,
Σ [dP _chem. (L) / dt] = k _react. · C _sub. · Σ [C _fusion (l)]
Based on chemiluminescence per unit time observed in the entire cell population: Σ [dP _chem. (L) / dt], the sum of the intracellular concentrations C _fusion (l) of the fusion protein in the entire cell population: Σ [C _fusion (l)] can be estimated.

  That is, a photoprotein (luciferase) gene is inserted as a “Tag”, in this case, a “reporter enzyme protein” gene so as to quantitatively reflect the expression level of a specific gene, in particular, the frequency of transcription. It is produced as a protein, and the enzyme activity of this “Tag” photoprotein (luciferase) is measured nondestructively. In particular, as a “reporter enzyme protein”, a luciferase gene encoding a photoprotein (luciferase) is inserted immediately after a sequence (IRES) capable of continuing transcription of a specific gene in the genomic DNA of the target plant. A transformant is obtained. About this transformant plant, each plant individual is self-crossed, and each T1 seed is harvested and selected as a strain that retains the above-described genetic recombination (line), in particular, “LucTag line” It is called.

  This “LucTag line” is a transformant plant, but the expression of a specific gene into which “LucTag” is inserted follows essentially the same control mechanism as that of the wild type plant. That is, the mechanism governing the expression induction, expression suppression, transcription, translation, protein folding, and degradation of the produced protein is the transforming plant of the “LucTag line”. It is essentially the same between wild-type plants. By utilizing this feature, while the plant of the “LucTag line” type transformant plant is grown, the temporal change in the expression level of a specific gene in the cell population of the plant body is traced. Information completely equivalent to the temporal change in the expression level of a specific gene in a cell assembly of a wild-type plant can be obtained.

The insertion of the “reporter enzyme protein” luciferase gene into the chromosomal DNA of a wild-type plant is a homologous chromosome pair (a · a) in which a specific gene is present in a “diploid” type chromosomal DNA, Occurs on either homologous chromosome. In other words, a transformed plant body that has been subjected to genetic recombination operation has a chromosome configuration of type (a ^* .a) consisting of a chromosome a ^* that has been genetically modified and a homologous chromosome a that has not been genetically modified. Is adopted. Therefore, when plant individuals (a ^* · a) of the obtained transformed plants are self-mated and T1 seeds of each plant individual are harvested, the T1 seeds obtained are (a ^* · a ^* ), (a ^* )・ A), (a ^* a ^* ), and (a ^* a) are selected from any of the four types of chromosome configurations. As the T1 seeds are lines of, the "LucTag line" type transformant plants, (a ^* · a ^*) with a chromosomal structure of the "homo-type", as well as, (a ^* · a) or (a -"Heterotype" having a chromosome structure of a ^* ) and seeds with different genotypes are mixed.

Even if the plant individual seeded with the seed of the “LucTag line” type transformant, germinated, and grown is, for example, a “homo type” having a chromosome configuration of (a ^* · a ^* ), self-mating In this case, allelic recombination partially occurs between homologous chromosome pairs, and thus there may be a difference in phenotype between individuals. That is, even with a “homotype” having a chromosomal configuration of (a ^* · a ^* ), the temporal change in the expression level of a specific gene to which “LucTag” is attached varies between plant individuals. The distribution is shown. Further, between the plants in the "homo-type" with chromosomal structure, and plants in the "heterozygous" with a chromosomal structure of the (a ^* · a) or (a · a ^*) of the (a ^* · a ^*) Among them, due to the difference in genotype, the temporal change in the expression level of a specific gene to which “LucTag” is attached is significant among plant individuals having different genotypes. May show differences.

  In addition, when the “reporter enzyme protein” luciferase gene is inserted into the chromosomal DNA of a wild-type plant, the insertion site may accidentally inhibit the translation of the C-terminal part of the coding region of the “specific gene”. In this case, the C-terminal deletion mutant protein lacking the C-terminal portion, not the wild-type protein originally produced by the “specific gene”, due to the expression of the inserted “specific gene”. Is produced. In some cases, the function of the wild-type protein is impaired or lowered in the C-terminal deletion mutant protein. For example, if the inserted “specific gene” is a self-regulating gene, its function is lost, so that “homotype” and “heterotype” There may be a difference in the temporal change of the expression level in the cell.

  When using the “LucTag line” type transformant plant, each genotype is examined after verifying the dispersion between individuals within the same genotype and the presence or absence of phenotypic differences between different genotypes. It is necessary to extract a common tendency (feature) of “temporal change in the expression level of a specific gene”.

  Specifically, first, a plurality of seeds of the “LucTag line” type transformant plant were sown, germinated and grown, each plant individual was measured using “LucTag”. Time series gene expression level data showing the “temporal change in the expression level of a specific gene” observed in the entire cell population of the plant individual is compared with each other, and the degree of dispersion among individuals within the same genotype As well as verifying the presence or absence of phenotypic differences between different genotypes. Next, when there is a significant difference in phenotype between different genotypes, a common “temporal change in expression level of a specific gene” within each genotype is extracted. It will be necessary.

  On the other hand, if there is a slight phenotypic difference between different genotypes, within the same “line”, the different “generic gene expression level” It is necessary to extract the tendency (feature) of “temporal change”.

  Using the above-mentioned “LucTag line” type transformant plant, in order to elucidate the tendency (characteristics) of the “temporal change in the expression level of a specific gene”, the above-described analysis work is required. Become. In order to perform the analysis work efficiently and while maintaining high validity, it is desired to develop the following numerical and statistical analysis methods.

That is, as described above, “time-series data (t _i , P (t _i )) indicating a quantitative variation (waveform) with time, which is highly likely to indicate some“ commonness / similarity ”. ”For time-series data (t _i) for the purpose of“ cluster analysis ”that classifies hierarchically according to the statistical height of“ commonality / similarity ”. , P (t _i )) ”as well as“ time series data (t) determined to have statistical “commonness / similarity” as a result of “cluster analysis”. _i , P (t _i )) ”for a subset (class), statistical data found between“ time-series data (t _i , P (t _i )) ”included in the subset (class) “Time-series” with statistical “averaging” that reflects common “similarity / similarity” Data _{(t i, P (t i} )) "statistical processing for obtaining a representative value, and the resulting" time-series data _{(t i, P (t i} )) "based on the representative value, the subset “A feature on the waveform” indicating a statistical “commonness / similarity” found between “time-series data (t _i , P (t _i ))” included in (class) “ Development of statistical analysis methods for the purpose of profiling is desired.

  The present invention solves the above-mentioned problems, and an object of the present invention is to extract features of time-series gene expression level data in an individual organism composed of multiple cells that is suitable for the above-mentioned “use”, and An object of the present invention is to provide an analysis method for the purpose of comparative classification based on characteristics of time-series gene expression level data between individual organisms, and an analysis apparatus based on the analysis method. In particular, an object of the present invention is to use the recombinant plant in an individual plant for a genetically modified transgenic plant that meets the above-mentioned “use” when using a “LucTag line” type transformed plant. To provide an analysis method for the purpose of comparative classification based on feature extraction of time-series gene expression data of genes, characteristics of time-series gene expression data between individuals, and an analysis apparatus based on the analysis method is there.

  The genomes of various organisms have been deciphered, and many experiments are currently underway to make the meaning (search for functions; annotation) of the information contained in the genomes. Among the experimental systems, mutants and cultured cells that can follow a specific gene expression in time series have been developed, and there is a movement to elucidate the function of genes using these. A number of mutants and cultured cells (lines) have been prepared, and research has been conducted to comprehensively explore gene functions.

  The time-series gene expression data obtained in such an experimental system draws complex waveforms unique to the line, and even in the same line, the waveform varies depending on the genotype. Therefore, in the conventional manual method, it is not easy to grasp the characteristics of the waveform of each line and compare the similarity between lines.

  The analysis method according to the present invention is developed as an analysis tool specialized in processing time-series gene expression data obtained in such an experimental system, grasping the characteristics of the lines, and performing comparative classification of each line. It has been done. From the obtained data, the waveform characteristics of each line are displayed in an easy-to-see manner, the genotype of the sample is automatically determined from the waveform characteristics, and a model waveform is created and analyzed for each. Finally, cluster analysis is performed using the model waveform of each line, and comparison classification of each line is performed.

  The analysis system according to the present invention has the following functions.

  The “measurement error correction” 14 receives numerical data obtained from an experimental measuring instrument, detects if there is a false positive in the value of gene expression, and corrects it.

  The “plurality waveform discrimination” 17 determines whether a plurality of different waveforms are mixed depending on the genotype or the like from the value of gene expression. This is because the waveform may differ greatly depending on the genotype, and it is necessary to distinguish between the two when processing the model waveform creation 7 described later.

  “Unify time axis” 18 unifies the time when the start time and measurement interval of the gene expression measurement experiment differ between lines.

  “Model waveform creation” 19 creates a waveform that most reflects the characteristics of the line from a plurality of samples of the line. The model waveform can be said to be standardized data by correcting genotyping errors and measurement time errors between lines.

  The “waveform decomposition” 20 assumes that a model waveform (complex curve) is a composite of known functions (simple curves such as a normal distribution), and decomposes them into a plurality of simple curves. By thinking in this way, complex waveforms, particularly peaks, can be seen more clearly. A profile of time, half width, and vertex is created for each separated curve.

  “Waveform comparison / classification” 23 performs cluster analysis based on the data standardized in the above waveform analysis 5 and statistically classifies similar lines.

That is, the method for analyzing time-series gene expression level data according to the present invention includes:
An analysis method for the purpose of comparative classification based on feature extraction of time series gene expression data in individuals derived from the same "line" and characteristics of time series gene expression data between individuals,
The time-series gene expression level data in the individual organism uses a luciferase gene as a `` reporter gene '' capable of non-destructively monitoring the expression level of the gene, and as a chronological change in `` luciferase enzyme activity '', Observed data,
The analysis method is
It is determined whether or not “measurement error” data caused by the measuring instrument is mixed in the time-series numerical data of “luciferase enzyme activity” obtained from the measuring instrument. A process of “measurement error correction”, which performs an operation of correcting the “measurement error” data to an “estimated value” based on the numerical data before and after that;
Comparing numerical data of time-series changes of “luciferase enzyme activity” for multiple organisms originating from the same “line”, hierarchical clustering is performed according to the level of similarity, and at least Performing “cluster analysis” operation to form one or two groups including the individual of “multiple types of waveform discrimination”;
Based on numerical data of time-series changes in the “luciferase enzyme activity” for biological individuals originating from the same “line”, observed in the individual at a desired elapsed time starting from the measurement start time It is converted into numerical data of time-series changes in predicted “luciferase enzyme activity”, and “time axis unification” is performed between the numerical data of time-series changes in “luciferase enzyme activity” in each individual. Process of “unification of time axis”;
Using the numerical data of the time-series changes in “luciferase enzyme activity” in each individual for which “time axis unification” has been made, a plurality of types determined to belong to the same group by the “plurality waveform discrimination” Create a “model waveform” that shows the time-series changes in the representative value of “luciferase enzyme activity”, reflecting the similarity of the time-series changes in “luciferase enzyme activity” in this group. The process of “model waveform creation”;
Based on a “model waveform” indicating a time-series change of a representative value of “luciferase enzyme activity” for a group consisting of a plurality of individuals determined to belong to the same group, a plurality of “unimodal waveform functions” are provided. Overlay, create a “synthetic waveform” that approximately shows the waveform characteristics of the “model waveform”, and for each of the “unimodal waveform function” values for the peak position, peak height, and half width A method of analyzing time-series gene expression level data, characterized by having a step of “waveform decomposition”.

In addition, the method for analyzing time-series gene expression level data according to the present invention includes:
Using the “model waveform” created in the process of “model waveform creation” for the biological group derived from each “line”,
Compare “model waveforms” of organisms derived from multiple “lines”, perform hierarchical clustering according to the level of similarity, and at least one group that includes multiple “lines”, or It is preferable that the method further includes a “waveform comparison / classification” step of performing two “cluster analysis” operations.

In the method for analyzing time-series gene expression level data according to the present invention,
The “unimodal waveform function” used in the “waveform decomposition” step is preferably a Lorentz function type waveform function.

In the process of “model waveform creation”,
In a group consisting of a plurality of individuals determined to belong to the same group, as the representative value of “luciferase enzyme activity” that reflects the similarity of time-series changes in “luciferase enzyme activity” in this individual group, the individual It is preferable to select the central value of the “luciferase enzyme activity” at each time in the group.

In the “Multiple Waveform Discrimination” process,
When comparing numerical data of time-series changes in the “luciferase enzyme activity” for multiple organisms that originate from the same “line”, each individual can be clustered according to the level of similarity. The similarity of the numerical data of the time-series changes in the “luciferase enzyme activity” between the two is based on the matrix matrix calculation method for calculating the “distance” between the numerical data of the time-series changes in the “luciferase enzyme activity”. It is desirable to use it.

that time,
As the inter-row matrix calculation method, a manhattan method is used,
It is preferable to use the ward method as the binding method for clustering based on the “distance” between the numerical data of the time-series changes in the “luciferase enzyme activity”.

In addition, the analysis system according to the present invention includes:
An analysis system that can be used for analysis for the purpose of comparative classification based on feature extraction of time-series gene expression data in organism individuals originating from the same "line" and characteristics of time-series gene expression data between individuals Because
The analysis system is a mechanism for performing analysis according to the method for analyzing time-series gene expression level data of the present invention.
It is determined whether or not “measurement error” data caused by the measuring instrument is mixed in the time-series numerical data of “luciferase enzyme activity” obtained from the measuring instrument. A “measurement error correction” mechanism that performs an operation to correct the “measurement error” data to an “estimated value” based on the numerical data before and after that;
Comparing numerical data of time-series changes of “luciferase enzyme activity” for multiple organisms originating from the same “line”, hierarchical clustering is performed according to the level of similarity, and at least A "multiple-type waveform discrimination" mechanism that performs a "cluster analysis" operation to form one or two groups, including
Based on numerical data of time-series changes in the “luciferase enzyme activity” for biological individuals originating from the same “line”, observed in the individual at a desired elapsed time starting from the measurement start time It is converted into numerical data of time-series changes in predicted “luciferase enzyme activity”, and “time axis unification” is performed between the numerical data of time-series changes in “luciferase enzyme activity” in each individual. "Unified time axis"mechanism;
Using the numerical data of the time-series changes in “luciferase enzyme activity” in each individual for which “time axis unification” has been made, a plurality of types determined to belong to the same group by the “plurality waveform discrimination” Create a “model waveform” that shows the time-series changes in the representative value of “luciferase enzyme activity”, reflecting the similarity of the time-series changes in “luciferase enzyme activity” in this group. The “model waveform creation” mechanism;
Based on a “model waveform” indicating a time-series change of a representative value of “luciferase enzyme activity” for a group consisting of a plurality of individuals determined to belong to the same group, a plurality of “unimodal waveform functions” are provided. Overlay, create a “composite waveform” that approximately shows the waveform characteristics of the “model waveform”, and for each of the “unimodal waveform functions”, the values of the peak position, peak height, and half width It is a time-series gene expression level data analysis system characterized by having a “waveform decomposition” mechanism.

In the time series gene expression level data analysis system according to the present invention,
Using the “model waveform” created in the process of “model waveform creation” for the biological group derived from each “line”,
Compare “model waveforms” of organisms derived from multiple “lines”, and perform hierarchical clustering according to the level of similarity, and at least one group that includes multiple “lines”, or The system configuration further includes a “waveform comparison / classification” mechanism that performs two “cluster analysis” operations.

The analysis system and analysis method according to the present invention particularly have the following effects.
・ When a large amount of time-series expression data is obtained, the lines can be compared and classified with quantitative data that has been objectively analyzed.
・ By performing waveform decomposition, a single waveform is viewed as a plurality of waveforms, so that the characteristics can be seen more remarkably.

  Hereinafter, the present invention will be described in more detail.

(1) “LucTag line” type transformed plant to be measured First, “temporal change in expression level of a specific gene” targeted by the analysis method of the present invention is a “diploid” type chromosomal DNA. Means that the expression level of a “specific gene” present as an allele in each homologous chromosome pair varies with time in a state in which the growth of the plant body is maintained. Specifically, it is not a “temporal change in the expression level of a specific gene” within a single cell of the plant body, but a part of the plant body, for example, a specific organ or tissue such as a leaf. Is a “temporal change in the expression level of a specific gene” in the entire cell population including a large number of cells.

The number of cells constituting the cell population N _total gradually changes with time in a state in which the growth of the plant body is maintained. For example, the cell population constituting an organ in which differentiation such as a leaf has been completed Then, the number of cells constituting the cell can be considered to be substantially constant in a time width of about several days. In organs such as leaves, the cells constituting this cell population are placed in the same environment as a whole. An organ such as a leaf, for example, a group of epidermal cells of a leaf, is a group of cells of the same type, but when viewed microscopically, there are slight differences depending on the site where the cells exist. Therefore, the number N _total of cells constituting the cell population can be expressed as the sum of the number N _subgr-i of several subsets, N _total = ΣN _subgr-i . Microscopically, the subset of cells is substantially the same type of cell group.

  At that time, microscopically, in a group of cells of substantially the same type, in consideration of “temporal change in the expression level of a specific gene” in the cell, the expression of the specific gene starts and the expression level The “temporal change” that reaches the maximum and eventually decreases can be considered to be substantially the same. On the other hand, even in a group of cells of substantially the same type that are placed in the same environment, the time at which the expression of the specific gene starts is not completely synchronized. It can be approximated according to a probability function that maximizes a certain time.

For example, like a transformed E. coli, a single clone, i.e., all cells, are genetically composed of completely the same type of cells, and in a culture medium, a recombinant gene that utilizes an inducer. Even when overexpression is performed, if the frequency of overexpression in the cell group is traced over time, the time of addition of the inducer is set to t = 0, and the number of cells that start overexpression gradually increases at a certain time. It reaches a maximum at t _max and then decreases. The time variation of the number Ns (t) of cells that start this overexpression is Poisson distribution: f (x) = k ^x · exp (−k) / x! It can be approximated as a probability distribution similar to. Once overexpressed, if the protein of the gene product is produced in the cell in a short period of time and then the protein is not degraded, the total amount of the protein in the entire cell group is the overexpressed cell. It is proportional to the sum of numbers. Total number of overexpressed cells: ∫Ns (t) dt is the number of cells Ns (t _i ) that start overexpression at each time t _i .
∫Ns (t) dt≈Σ1 / 2 · {Ns (t _i ) + Ns (t _{i + 1} )} · (t _{i + 1} −t _i )
And can be expressed approximately. The time change of the number N (t) of cells in which overexpression starts is determined according to the Poisson distribution. Ns (t) = Na × {t _max ^t · exp (−t _max ) / t! }
If so denoted, overexpressed cell number of total: ∫Ns (t) dt is a t = t _max, the increase rate Ns (t) represents the maximum, is reached t = 2t _max, the saturation tendency It becomes a monotonically increasing function. Correspondingly, considering the total amount P _pro (t) of the protein in the entire cell group, the rate of increase at each time: dP _pro (t) / dt is
dP _pro (t) / dt∝N (t)
dP _pro (t) / dt = P _pro a × {t _max ^t · exp (−t _max ) / t! }
When t = t _max , the rate of increase: dP _pro (t) / dt shows a maximum, and when t = 2t _max is reached, a monotonically increasing function showing a saturation tendency.

Once over-expressed, the gene product protein is produced in the cell in a short period of time, and when the protein is subsequently enzymatically degraded, the apparent rate constant k of the enzymatic degradation reaction of the protein. _{According to dig.} , the concentration of the protein in the cell decreases. At that time, the number Ns (t) of cells initiating overexpression is determined according to the Poisson distribution.
Ns (t) = Na × {t _max ^t · exp (−t _max ) / t! }
Can be expressed as 成分 P _{pro . +} (T) / ∂t = P _pro a × {t _max ^t · exp (−t _max ) / t! }
It becomes. Correspondingly, the total amount P _pro (t) of the protein in the whole cell group is
P _pro (t) = ∫ [(∂P _{pro. +} (S) / ∂s) · exp {−k _dig. (Ts)}] · ds
≈Σ1 / 2 · [{(∂P _{pro. +} (T _i ) / ∂t) · exp {−k _dig. (T−t _i )} + {(∂P _{pro. +} (T _{i + 1} ) / ∂t) · exp {−k _dig. (T−t _{i + 1} )}] · (t _{i + 1} −t _i )
And can be expressed approximately. At that time, the total amount P _pro (t) of the protein in the whole cell group can be secondarily differentiated so as to show a maximum between the time t _max and the time t _dig. ≡ {1 / k _dig. } _. It can be approximated as a continuous function.

Regarding “a temporal change in the expression level of a specific gene” in a whole cell population including a large number of cells constituting a specific organ or tissue such as a part of a plant, for example, a leaf, When the number N _{total of} cells constituting a cell population can be expressed as the sum of the number N _subgr-i of several subsets of cells of the same kind, N _total = ΣN _subgr-i , the number N _{subgr- The} “temporal change in the expression level of a specific gene” in the same type of cell group _i can be approximated by the above-described model.

That is, in the same cell group having the cell number N _subgr-i , the cell number Ns _i (t _i ) at which expression of the “specific gene” starts at time t _i is determined according to the Poisson distribution.
Ns _i (t) = Na _i × {t _max-i ^t · exp (−t _max-i ) / t! }
When, approximately can notation, increase component of amount of protein from cells to initiate this _{expression, ∂P pro + -. I (} t) / ∂t = P pro-i a × {t max-i t Exp (-t _max-i ) / t! }
It becomes. Correspondingly, the total amount P _pro-i (t) of the protein in the entire cell group of this subset of cell numbers N _subgr-i is
P _pro-i (t) = ∫ [(∂P _{pro. + -I} (s) / ∂s) · exp {−k _dig. (Ts)}] · ds
≈Σ1 / 2 · [{(∂P _{pro. + -I} (t _i ) / ∂t) · exp {−k _dig. (T−t _i )} + {(∂P _{pro. + -I} (t _{i +1} ) / ∂t) · exp {−k _dig. (T−t _{i + 1} )}] · (t _{i + 1} −t _i )
And can be expressed approximately. In this case, the total amount P _pro-i (t) of the protein in the entire cell group constituting the subset is between the time t _max-i and the time t _dig. ≡ {1 / k _dig. } _. It can be approximated as a continuous differentiable function that exhibits a maximum.

The target cell population is a collection of several subsets of the same type of cells, and the total amount P _pro (t) of the protein expressed from the “specific gene” in the entire cell population is: P _pro (t) = ΣP _pro-i (t)
Can be expressed as: In this case, in the temporal change of the total amount P _pro (t) of the protein expressed from the “specific gene” in the entire cell population, the portion of the cell group constituting a plurality of subsets in a certain time zone The total amount P _pro-i (t) of the proteins in the assembly shows a maximum and may overlap each other.

On the other hand, the total amount P _pro-i (t) of the protein in the entire cell group constituting the subset is a _maximum between the time t _max-i and the time t _dig. ≡ {1 / k _dig. } _. As shown in FIG. 4, at least the apparent full width at half maximum Δt _hw of the peak is extremely large compared to t _max-i or t _dig. ≡ {1 / k _dig. } _. It can never be smaller. Therefore, even in the temporal change of the total amount P _pro (t) of the protein expressed from the “specific gene” in the entire cell population, the apparent full width at half maximum Δt _hw of the peak is t _max-i or t _{dig. . ≡ {1 / k dig.} } as compared to, is impossible to become extremely small. When observing the total amount of protein expressed from a “specific gene” in the entire cell population of interest, compared to the t _max-i or t _dig. ≡ {1 / k _dig. }, When an “extremely sharp peak” showing a full width at half maximum Δt _hw that is extremely small is observed, this “extremely sharp peak” is likely to be a “noise peak” caused by the measurement system. Extremely expensive.

  Certainly, the possibility that many cells do highly synchronous “expression” cannot be ruled out, but this type of highly synchronous “expression” is not limited to certain “such as a“ heat-shock ”protein. Limited to genes encoding very few types of maintenance proteins that need to respond instantaneously to "stress".

In many cells, when the “cell cycle” proceeds with extremely high synchrony (synchronous division occurs), each stage of the “cell cycle” (S phase, G ₂ phase, M phase) The gene product (protein) of the gene expressed only in the (phase) shows a change in concentration according to the period of the “cell cycle”. In consideration of the fact that cell populations constituting plant bodies, particularly cell populations constituting already differentiated organs and tissues, generally have a low mitotic index and range from several to 20%. Then, the possibility of observing a “steep peak” due to the above-mentioned synchronized splitting cannot be excluded, but is extremely rare.

  In the present invention, as a means for detecting “temporal change in the expression level of a specific gene”, a form utilizing a “LucTag line” type transformant plant is selected.

  Specifically, for a specific gene present on the genomic DNA of a plant, a ribosome containing a region encoding a protein amino acid sequence encoded by the gene is continued so as to continue transcription. A region (ORF portion) encoding the photoprotein (luciferase) of the luciferase gene is inserted immediately after the sequence (IRES) that enables binding to and linked. This chimeric gene has a promoter sequence that causes transcription of a specific gene in the upstream portion, and transcription is similarly initiated by an endogenous factor that induces expression of the specific gene. In that case, the obtained mRNA becomes translatable into a peptide chain of a photoprotein (luciferase) encoded by the inserted luciferase gene. That is, translation is performed based on the mRNA, and when folding proceeds into a protein, an active photoprotein (luciferase) is produced in the cell as a reporter protein. The expression level corresponds to the expression level of the target protein encoded by the target specific gene in the wild-type plant.

  In a specific example to be described later, Arabidopsis thaliana is created using a technique of inserting a DNA fragment having a known base sequence containing a reporter gene into the chromosomal DNA of Arabidopsis thaliana using the IRES. The transformant is used. At this time, a condition is selected in which insertion of a known DNA fragment occurs randomly in any one of a plurality of sites having a base sequence of the IRES present in chromosomal DNA. Accordingly, different types of transformants can be obtained depending on which of the plurality of sites has a known DNA fragment inserted therein.

  In particular, as a reporter gene, a luciferase gene encoding a firefly-derived luciferase, which has been conventionally verified to produce a protein having an enzyme activity when recombinantly expressed in a plant body, is used. In addition, firefly-derived luciferase is a recombinantly expressed protein having enzyme activity in cells constituting various organs of plants, for example, roots, stems, and leaves, as verified in firefly grass, for example. Produced. When the substrate of the luciferase derived from the firefly is absorbed from the root of the plant body and supplied to the cells of each organ and tissue of the plant body via the conduit, the enzymatic activity of the recombinantly expressed luciferase in the cell, The substrate luciferin is converted to oxyluciferin, and blue chemiluminescence derived from the oxyluciferin occurs in an equal amount.

  A plurality of types of transformants having different insertion sites are self-mated to harvest T1 seeds of each plant individual. This T1 seed is sown, and in the grown plant body, a line in which recombinant expression of firefly-derived luciferase is confirmed is self-mated and T2 seed is collected. In a specific example, T2 seeds collected from a line showing the phenotype of the recombinant expression of firefly-derived luciferase are sown, and the grown plant body is measured.

  In addition, the method disclosed in the following literature is applied to the method of creating the “LucTag line” type transformant of Arabidopsis thaliana.

References:
The Plant Journal (2003) 35, 273-283
Gene trapping of the Arabidopsis genome with a firefly luciferase reporter
Yoshiharu Y. Yamamoto, Yumi Tsuhara, Kazuhito Gohda, Kumiko Suzuki and Minami Matsui

  In the “LucTag line” type transformant plant created by the above method, since the inserted DNA fragment itself does not have a promoter sequence, the expression of a gene that induces endogenous expression in the plant. Along with this, transcription of the “reporter gene” inserted downstream is carried out. Furthermore, when two or more “reporter genes” are inserted into the chromosomal DNA, the two “reporter genes” may be transcribed. Since the amount of protein of recombinantly expressed luciferase accompanying the transcription of the luciferase gene of the “reporter gene” reflects the expression level of one specific gene, The transformant in which the “gene” is inserted is selected as the “LucTag line” type transformant used in the present invention.

That is, at the stage of line formation, the chromosomal DNA of the plant growing from the sowed T1 seed is collected, the chromosome a ^{* in} which the “reporter gene” has been recombined, and the homologous chromosome in which the gene recombination has not been performed. It is confirmed that it is a “heterotype” having a chromosome structure of (a ^* · a) type consisting of a. In addition to confirming the phenotype, the “LucTag line” -type transformed plant that has also been confirmed to be “hetero-type” is self-crossed to collect T2 seeds. T2 seeds collected from this “heterotype” transformed plant have a “homotype” having a chromosome configuration of (a ^* · a ^* ), and (a ^* · a) or (a · a ^* ) A “heterotype” having a chromosomal configuration and a “wild type” having a chromosomal configuration of (a · a) are mixed with seeds having different genotypes.

When self-mating is carried out, genetic recombination occurs between homologous chromosomes, resulting in crossover between gene linkage groups. In T2 seeds, (a ^* · a) or (a · a ^* ) Even in “heterotype” having a chromosomal structure, there is a difference in “genotype”. At this time, there may be a difference in expression frequency due to the structure of the gene linkage group on the chromosome where the “target gene” in which the target “reporter gene” has been recombined is present. Alternatively, when expression occurs in parallel in the “target gene” and its allele, there may be a difference in the frequency ratio of the expression due to the configuration of the gene linkage group.

To explain this point more specifically, when the composition of the gene linkage group of the chromosome of the T1 seed is (A ₁ A ₂ ^* , a ₁ a ₂ ), At least two types of (A ₁ A ₂ ^* , a ₁ a ₂ ) or (A ₁ a ₂ , a ₁ A ₂ ^* ) are mixed. At that time, when the “target gene” A ₂ ^* and its allele a ₂ are expressed in parallel, in the configuration of (A ₁ A ₂ ^* , a ₁ a ₂ ), the expression frequency is “ “Target gene” A ₂ ^* > allele a ₂ , but in the configuration of (A ₁ a ₂ , a ₁ A ₂ ^* ), the expression frequency is “target gene” A ₂ ^* <allele a ₂ A phenomenon is also conceivable.

Similarly, in the “homotype” T2 seed, the composition of the gene linkage group of the chromosome is at least (A ₁ A ₂ ^* , A ₁ A ₂ ^* ) or (A ₁ A ₂ ^* , a ₁ A _2. ^* ) The two types are mixed. Even between these two types, the expression frequency of the “target gene” A ₂ ^* may be slightly different.

  Therefore, even in the “LucTag line” based on one T1 seed, in particular, among the plant bodies grown from a plurality of T2 seeds, there are specific types of “homotype” and “heterotype” genotypes. There may be a difference in the temporal change in the amount of the “reporter gene” product, recombinantly expressed luciferase protein associated with gene expression. Furthermore, there is a possibility that a microscopic difference reflecting the “genotype” may be found among plant groups roughly classified into “homotype” and “heterotype”.

When the structure of the chromosomal gene linkage group in the wild type plant body used for the creation of the “LucTag line” is (A ₁ A ₂ , a ₁ a ₂ ), the insertion of the luciferase gene of the “reporter gene” When the allele “A2” is present on the homologous chromosome, the transformed plant has a genotype of (A ₁ A ₂ ^* , a ₁ a ₂ ). T1 seeds collected by self-mating this transformed plant are “homotype” (A ₁ A ₂ ^* , A ₁ A ₂ ^* ), (A ₁ A ₂ ^* , a ₁ A ₂ ^*). ), (A ₁ A ₂ ^* , a ₁ A ₂ ^* ); “Hetero” (A ₁ A ₂ ^* , a ₁ a ₂ ), (A ₁ A ₂ ^* , A ₁ a ₂ ), (a _{1 ^{a 2 *, a 1 a 2}} ), (a 1 a 2 *, a 1 a 2); "wild-type" of _{(a 1 a 2, a 1} a 2), (a 1 a 2, a 1 a 2 ), (A ₁ a ₂ , a ₁ a ₂ ). Thereafter, in the plant grown from the T1 seed, the product of the “reporter gene” and those capable of producing the recombinantly expressed luciferase protein are the “homotype” (A ₁ A ₂ ^* , A _{1 ^{a 2 *), (a 1}} a 2 *, * a 1 a 2), (a 1 a 2 *, a 1 a 2 *); "heterozygous" in ^{_{(a 1 a 2 *, a}} 1 a 2) , (A ₁ A ₂ ^* , A ₁ a ₂ ), (a ₁ A ₂ ^* , A ₁ a ₂ ), (a ₁ A ₂ ^* , a ₁ a ₂ ).

Next, when the composition of the gene linkage group of the chromosome of the T1 seed is “heterotype” (A ₁ A ₂ ^* , a ₁ a ₂ ), in the plant “LucTag line” grown from the T1 seed In T2 seeds collected by self-mating, “homotype” (A ₁ A ₂ ^* , A ₁ A ₂ ^* ), (A ₁ A ₂ ^* , a ₁ A ₂ ^* ), (a ₁ A ₂ ^* , a ₁ A ₂ ^* ); “heterotype” (A ₁ A ₂ ^* , a ₁ a ₂ ), (A ₁ A ₂ ^* , A ₁ a ₂ ), (a ₁ A ₂ ^* , A ₁ a ₂ ), (a ₁ A ₂ ^* , a ₁ a ₂ ); “wild type” (A ₁ a ₂ , A ₁ a ₂ ), (A ₁ a ₂ , a ₁ a ₂ ), (a ₁ a ₂ and a ₁ a ₂ ). Similarly, when the composition of the gene linkage group of the chromosome of the T1 seed is “heterotype” (a ₁ A ₂ ^* , A ₁ a ₂ ), in the plant body “LucTag line” grown from the T1 seed, In T2 seeds collected by self-mating, “Homo” (A ₁ A ₂ ^* , A ₁ A ₂ ^* ), (A ₁ A ₂ ^* , a ₁ A ₂ ^* ), (a ₁ A ₂ ^* , A ₁ A ₂ ^* ); “Hetero” (A ₁ A ₂ ^* , a ₁ a ₂ ), (A ₁ A ₂ ^* , A ₁ a ₂ ), (a ₁ A ₂ ^* , A ₁ a ₂ ), (a ₁ A ₂ ^* , a ₁ a ₂ ); “wild type” (A ₁ a ₂ , A ₁ a ₂ ), (A ₁ a ₂ , a ₁ a ₂ ), (a ₁ a ₂ , A ₁ a ₂ ).

On the other hand, when the structure of the gene linkage group of the chromosome of T1 seed is “heterotype” (A ₁ A ₂ ^* , A ₁ a ₂ ), In T2 seeds collected by self-mating, “homotype” (A ₁ A ₂ ^* , A ₁ A ₂ ^* ); “heterotype” (A ₁ A ₂ ^* , A ₁ a ₂ ); It may include a combination of “wild type” (A ₁ a ₂ , A ₁ a ₂ ). Further, when the composition of the gene linkage group of the chromosome of T1 seed is “heterotype” (A ₁ A ₂ ^* , a ₁ A ₂ ), in the plant body “LucTag line” grown from the T1 seed, In T2 seeds collected by self-mating, “homo-type” (A ₁ A ₂ ^* , A ₁ A ₂ ^* ); “hetero-type” (A ₁ A ₂ ^* , a ₁ A ₂ ); It may contain a combination of “wild type” (a ₁ A ₂ , a ₁ A ₂ ).

  Thus, even if it is roughly classified as “heterotype”, the composition of the collected T2 seeds has a large difference due to the difference in the composition of the gene linkage group of the chromosome of the T1 seed to be lined. It will be a thing. In verifying the presence or absence of differences between the “LucTag lines”, which are broadly classified as “heterotypes”, the number of plants grown from T2 seeds for measurement exceeds the number of combinations described above for each line. It is necessary to be.

Taking into account that point, in the specific examples below, plants grown from the T1 seeds in "LucTag line", by self-pollination, of the T2 seeds are collected, at least, 2 ⁴ (16) As a set, plants that have been sown and grown are associated with the expression of a specific gene, the product of a “reporter gene”, the presence or absence of a recombinantly expressed luciferase protein, and the recombinantly expressed luciferase protein. The amount of time change is measured.

In each line, contained in the plants grown from T2 seed, when the "wild-type" population is large beyond expectation, it is further possible to add 2 ⁴ (16) one set preferable. That is, even under the same experimental conditions, when the expression of a specific gene in each plant individual follows a “Poisson distribution” type random variable as described above, Since there is a possibility that the peak height may vary, the number of individuals producing the product of the “reporter gene” and the recombinantly expressed luciferase needs to be a certain number or more.

  In each line to be measured, in each plant grown from the T2 seed, an identifier of the “line” name given to the T1 seed that is the origin, and an identifier given to the T2 seed An ID for identifying an individual of “line” name- “branch number” type is attached. Further, the identifier of the “line” name attached to the T1 seed is attached to the portion where “LucTag” is inserted, that is, the identifier that specifies the “specific gene” and the T1 seed. The identifier is combined with the “specific gene” name—the “branch number” type.

  Each T1 seed has been confirmed to be recombinantly expressed at the stage where it has been sown and grown to the presence of a target “reporter gene” product and a recombinantly expressed luciferase protein. A “phenotype” of the presence or absence of a luciferase protein has been identified. Of course, at this stage, the “line” in which the presence of recombinantly expressed luciferase protein is not confirmed in the plant grown from the T1 seed is determined not to be a “LucTag line” type transformed plant.

(2) for measuring object "LucTag line" type transformation measured temporal variation of the "luciferase enzyme activity" in a plant individual "LucTag line" type transformation each line of plants, 2 ^four the T2 seeds (16 ) As a set, seeds are seeded one by one in each hole of a 96-well plate, and after germination is confirmed, “luciferase enzyme activity” is measured for the seedlings at predetermined time intervals. Specifically, the substrate luciferin germinated on an aqueous medium containing a predetermined concentration of the substrate luciferin and absorbed from its root is expressed recombinantly in the cell population constituting the above-ground organ of the seedling body. The intensity of blue chemiluminescence derived from enzymatically converted oxyluciferin is measured at predetermined time intervals based on the enzymatic activity of the luciferase.

A spectrum of blue chemiluminescence derived from oxyluciferin is known, and the light intensity of chemiluminescence having a narrow wavelength width including the peak wavelength of the chemiluminescence is measured. The amount of recombinantly expressed luciferase protein is proportional to the amount of expression of the “specific gene” of interest, so that the luciferase that is recombinantly expressed in the time zone when the expression of the “specific gene” is at a low level. The amount of protein remains at a low level. That is, in that time zone, the light intensity of blue chemiluminescence derived from the observed oxyluciferin is also low. The photon counting method is used to achieve sufficient measurement sensitivity even in the measurement of the light intensity of such extremely weak chemiluminescence. In the “photon counting method”, the photons received by the measuring device are counted one by one, and the time required from when the measurement is disclosed until a predetermined number of photons are counted is measured. The predetermined number of photons, that is, the time to reach the threshold photon number: N _{p-th .} : Δt _obs. Is used as the received light intensity, ie, the number of photons received by the measuring device per unit time: dN _photon / dt And
N _p-th. ≈ {dN _photon / dt} × Δt _obs.
Can be approximated. Therefore, the received light intensity: dN _photon / dt is
dN _photon / dt≈N _p-th. / Δt _obs.
Is expressed as:

However, when the received light intensity, that is, the number of photons received by the measuring device per unit time: dN _photon / dt is substantially “0”, the threshold number of photons reaches N _p-th. Time: Δt _{obs. Becomes} “∽”. Therefore, in actuality, when the number of photons N _p-gate received by the measuring device during the predetermined time width: Δt _gate does not reach the above-mentioned threshold number of photons: N _p-th. Light intensity received: dN _photon / dt is
Approximate as dN _photon / dt≈N _p-gate / Δt _gate .

The actual “photon counting” type measurement system uses a preset very short time width: δt _k and the number of photons received during that time: N _p-obs. (Δt _k ) _. , the integrated _{value:. ΣN p-obs (δt} k) is a threshold number of _photons:. N _p-th determines whether more than, at the time of exceeding, and a Δt _obs = Σδt _k.. The measurement sensitivity of the light receiver itself, the very short duration: when receiving several photons .DELTA.t _k, so as to be counted it with high accuracy, is set to a very high sensitivity. Specifically, the very short duration: in .DELTA.t _k, the threshold number of _photons:. If N _p-th 1/10 approximately of the light incident exceeding the number of photons is made greater than "measurable limit" In many cases, the measurement sensitivity is set so as to be

Temporarily, while measuring the light intensity of chemiluminescence, it is a very short time width, but compared with the light intensity of chemiluminescence that should be observed originally, it is remarkably caused by "pulsed stray light". When a large number of photons are incident on the light receiver, the extremely short time width: δt _k exceeds the “measurable upper limit”. Once the receiver reaches a state where the "measurable upper limit" is exceeded, the measurement is stopped thereafter, and a mode in which temporary damage to the photosensitive surface due to "light incidence exceeding the limit" is recovered is entered. . At that time, although the actual measurement is not completed, the “provisional measurement result” indicates the light intensity value of the “measurable upper limit” set as the measurement sensitivity. In fact, "pulsed stray light" is incident, very short duration: during .DELTA.t _k, the number of photons "Measurable limit": at the time of exceeding the _{N p-obs.LIMIT,} actually the value of _{N p-obs.current (δt k)} / δt k,; divided by .DELTA.t _k: counted once was the number of photons: the _{N p-obs.current (δt k)} , the very short duration Output as “provisional measurement result”. Of course, "tentative measurement result", the light intensity of "measurable limit": has a value greater than _{N p-obs.LIMIT} / δt _k.

When performing an analysis based on the present invention, after removing in advance the “wrong measurement result” caused by the “measurement error” of the receiver system, such as the above-mentioned “pulsed stray light” incidence, At that time (t _i ), the “measurement error corrected data” indicating the temporal change of the “luciferase enzyme activity”, which is estimated by “estimated” “measurement result” that is supposed to be observed and supplemented with the estimated value, is used. .

In addition, for each plant body (sample): S _Plant (m) to be measured, “original data” obtained by measuring the temporal change in “luciferase enzyme activity” is actually at a certain time (t _0-m ). Measurement is started, and then blue chemiluminescence light derived from oxyluciferin, which is an index indicating “luciferase enzyme activity” at the target time interval: Δt _interval in that order (t _im ). Measured value of intensity: P _lum.obs.-m (t _im ) measured in time series. That is, “data” in which a set of (t _im , P _lum.obs.-m (t _im )) of the time when actual measurement is performed and the measured value of the light intensity of chemiluminescence is recorded in time series. It has a format. The time (t _im ) when this actual measurement was made is in the format of HH: MM: SS [AM / PM].

(3) Processing of measurement data of temporal change of “luciferase enzyme activity” of each plant body (3-1) Conversion of time display to elapsed time display (“elapsed time conversion” processing)
Each plant (sample) to be measured: S _Plant (m) germinates after sowing, and the temporal change in “luciferase enzyme activity” is measured under predetermined growth conditions. For each "LucTag line", the T2 seeds as a ^two-four (16) a pair, are aligned with the time to start the experiment. Based on the start time (t _start-m ) of this experiment, the elapsed time up to the time (t _im ) when the actual measurement was made: t _m [i] = (t _im −t _start-m ) The time display type time-series measurement data (t _im , P _lum.obs.-m (t _im )) is changed from the elapsed time display type time series data (t _m [i], P _lum.obs.-m (T _m [i])).

(3-2) “Measurement error correction” process that eliminates “incorrect measurement results” caused by “measurement errors” in the receiver system and supplements with “estimated values” As described above, “luciferase enzyme activity ”Is predicted to be expressed by a continuous function that can be second-order differentiated: P _lum.-m (t _m ). In particular, the peak does not show a “pulse-like” sudden increase or decrease.

  On the other hand, the “wrong measurement result” due to the “measurement error” of the receiver system has a value that is at or exceeds the “measurable upper limit” of the receiver system. Compared with the “measurement result with no error” in the measurement time, the “pulse-like” increases and decreases abruptly.

As shown in the graph [A2] illustrated in FIG. 4, the “wrong measurement result” caused by the “measurement error” of the receiver system is a “pulse-like” often seen in so-called “spike noise”. Give a peak. Each measurement time interval: Δt _interval is at least sufficiently greater than the inactivation of recombinantly expressed “luciferase” or the reciprocal of the rate constant of the process undergoing degradation: t _dig. ≡ {1 / k _dig. } When the time constant t _dig. ≡ {1 / k _dig. }, The “pulse-like” peak generated with a time constant much shorter than the time constant t _dig. It can be determined as an “incorrect measurement result” due to an “error”.

“Luciferase enzyme activity” data at a certain time t _m [i]: _Is P _lum.obs.-m (t _m [i]) an “incorrect measurement result” due to a “measurement error” of the receiver system? The following conditions can be used as criteria for determining whether or not.

Specifically, the light intensity of the above “measurable upper limit”: P _lum-LIMIT ( _{≡N p-obs.LIMIT} / δt _k ) is used as a reference, and the light intensity is 1/10 or less. In a situation where the intensity is changing, the peak indicating a sudden increase / decrease of “pulse” is selected and removed as a measurement error.

Condition (1-0): The tip of the “peak” exists at time t _m [i].

{P _lum.obs.-m (t _m [i]) − P _lum.obs.-m (t _m [i−1])} × {P _lum.obs.-m (t _m [i + 1]) − P _lum.obs.-m (t _m [i])} <0
Condition (1-1): It is estimated that the light intensity is gradually changing.

| P _lum.obs.-m (t _m [i + 1])-P _lum.obs.-m (t _m [i-1]) | <{(P _lum.obs.-m (t _m [i + 1]) ) ^1/2 + (P _lum.obs.-m (t _m [i-1])) ^1/2 }
In the graph [A2] illustrated in FIG.
| P _lum.obs.-m (t _m [i + 1])-P _lum.obs.-m (t _m [i-1]) | <20
As a condition.

  Condition (1-2): In a situation where the light intensity is slowly changing, it is estimated that the range of the “statistically acceptable” dispersion is exceeded.

| P _lum.obs.-m (t _m [i])-1/2 · {P _lum.obs.-m (t _m [i + 1]) + P _lum.obs.-m (t _m [i-1] )} | <1/2 · {P _lum.obs.-m (t _m [i + 1]) + P _lum.obs.-m (t _m [i-1])}
_{P lum.obs.-m (t m [} i]) was measured and the result is, when following the Poisson distribution with a mean value _{{P lum.-m (t m} [i])}, the average value _{P lum.- Assuming that _m (t _m [i])} is an average value of _previous and subsequent measurement results, an allowable difference between the actual measurement value and the average value is the average value {P _lum.-m (t _m [i] )}

In the graph [A2] illustrated in FIG.
| P _lum.obs.-m (t _m [i])-1/2 · {P _lum.obs.-m (t _m [i + 1]) + P _lum.obs.-m (t _m [i-1] )} | <| Div. _max |
At that time, | Div. Select _max | = 100;
As a condition.

When the above conditions (1-0) to (1-2) are satisfied, a measurement error “peak” existing at the time t _m [i]: P _lum.obs.-m (t _m [i]) Instead of the value of 1/2 · {P _lum.obs.-m (t _m [i + 1]) + P _lum.obs.-m (t _m [i-1])} " _Completed time series data (t _m [i + 1], P _{lum.corrected.-m} (t _m [i])).

  When this “measurement error correction” is performed, the “measurement error correction” and “time series data” obtained are caused by the “measurement error” of the optical receiver system as shown in the graph [A3] illustrated in FIG. The “incorrect measurement result” is removed, and the waveform of “luciferase enzyme activity” that changes slowly and temporally is superimposed with the normal and statistically acceptable “measurement variation” It becomes.

(4) Determination of the presence or absence of expression of “luciferase gene” of “reporter gene” based on temporal change of “luciferase enzyme activity” of each plant body,
As shown in graph [A3] illustrated in FIG. 5, in the target “line”, the “luciferase enzyme activity” observed in the entire cell population of the measurement part of each plant is rising and falling. It is presumed to show a gradual change requiring 2 hours or more.

  In consideration of this point, “measurement error correction” indicates measurement results every two hours from time-series data that has been measured every hour as shown in the graph [A3] and has been subjected to “measurement error correction”. "Create time series data (sub set).

  That is, the time series data [A3] that has been subjected to the “measurement error correction” shown in FIG. 2 or the time series data [B2] that has been subjected to the “measurement error correction” process shown in FIG. Based on this, the presence or absence of expression of the “reporter gene” “luciferase gene” is determined.

  In addition to the above-mentioned “pulse stray light”, the measurement system also receives a “background / ground” type “stray light component” which is weakly incident but is continuously incident. Therefore, it is determined whether or not a clear light intensity level reflecting the temporal change in chemiluminescence, which indicates “luciferase enzyme activity”, has been continued for a certain period of time, and the “luciferase gene” of the “reporter gene” is determined. It is determined whether or not “

Specifically, since the “luciferase enzyme activity” observed in the entire cell population of the measurement part of each plant body is estimated to show a gradual change requiring 2 hours or more for both rising and falling, In the time-series data (sub-set) after “measurement error correction” indicating the measurement results every two hours, the light intensity at three times: P _{lum.corrected.-m} (t _m [i] ) Exceeds a certain level, it is determined that the “luciferase gene” of the “reporter gene” is expressed. The criterion that “significantly the light intensity due to chemiluminescence is measured” is
In time-series data [B2] that has been subjected to “measurement error correction” processing,
Time satisfying P _{lum.corrected.-m} (t _m [i]) ≧ 30: It is selected that t _m [i] continuously _extends over three time points, ie, 4 hours.

  The intensity of the “background” type “stray light component” that is continuously incident is observed in the “wild type” plant in which the “reporter gene” “luciferase gene” is not introduced under the same conditions. It is the same level as the light intensity. It is preferable to set the lower limit of the significance level to at least 2 to 3 times the intensity level of the “background” type “stray light component”.

Plants grown from T2 seeds that do not satisfy the above-mentioned criteria of “significantly, the light intensity caused by chemiluminescence is measured” have the “line” represented by (a ^* · a) or ( When it is derived from a T1 seed having a “heterotype” genotype having a chromosome configuration of a · a ^* ), it is normally considered to be a “wild type” having a chromosome configuration of (a, a).

Actually, the whole cell population of the measurement part of each plant body can be expressed as the sum of the number N _subgr-i of several subsets, N _total = ΣN _subgr-i, and the total number of subsets is In many cases, the cell number N _{subgr-i of} this subset is small, and when the temporal timing of its expression is completely dispersed, the light intensity measured over a long period of time: P _{lum.corrected.} The numerical integration value of _-m (t _m [i]) is
∫P _{lum.corrected.-m} (t) dt
≈Σ1 / 2 · {P _{lum.corrected.-m} (t _m [i + 1]) + P _{lum.corrected.-m} (t _m [i])} · (t _m [i + 1] −t _m [i])
As a whole, the integrated value exceeds a certain level, but the individual measurement time (t _m [i]) may not satisfy the above criteria.

Specifically, the numerical values of the measured light intensity: P _{lum.corrected.-m} (t _m [i]) at four consecutive measurement times are substantially as 29, 45, 45, and 29. In the case where it is 30 or more, and such four measurement time zones are dispersed over a long period of time, the above criteria are not “formally” satisfied. It is not considered that “the light intensity due to chemiluminescence is measured significantly”. That is, in the plant body, it is determined that there is no expression of the “specific gene” in which the “luciferase gene” of the “reporter gene” is inserted.

On the other hand, the numerical value of the measured light intensity: P _{lum.corrected.-m} (t _m [i]) at a single point and three consecutive measurement times over a long time is 31, 30, 31 It is 28, 31, 30, 31, and 27 in five consecutive measurement times including before and after that, and it is actually a case where it is unknown whether or not it will be 30 or more. However, if the above criteria are “formally” satisfied, it is considered that “the light intensity caused by chemiluminescence is significantly measured”. That is, it is determined that the “specific gene” into which the “luciferase gene” of the “reporter gene” is inserted is expressed in the plant body.

  The determination result (phenotype) of the presence or absence of “gene expression”, which is determined by the above-mentioned criteria for each measured plant body (sample) for the target “line”, is stored in the line table 31.

In case classification (branching) based on “presence / absence of expression”, for each target “line”, each plant body (sample) to be measured: S _Plant (m) has 2 ⁴ T2 seeds ( 16) As a set, for all of them, the determination result (phenotype) of presence or absence of “gene expression” is obtained, and at least one plant body (sample): S _Plant (m) is “gene expression” If it is determined that there is, analysis operations of “model waveform formation” and “line feature grasp” are performed.

Conversely, the target for each "line", the plants are measured (sample): S _Plant (m) is then the T2 seeds 2 ⁴ (16) and a set, for all, The determination result (phenotype) of presence or absence of “gene expression” is obtained, and if any plant body (sample): S _Plant (m) is determined to have no “gene expression”, “model waveform formation”, “ Do not perform the analysis operation of “Characteristics of line”.

  Accordingly, in the target “line”, each plant (sample) to be measured is a plant that does not clearly show a significant peak in “luciferase enzyme activity” as the whole cell population of the measurement part. In some cases, it is difficult to identify the peak and extract the features, but using this criterion, it is possible to eliminate those that lack this kind of “characteristic”.

By setting this standard, when the whole cell population of the measurement part of each plant body can be expressed as the sum of the number N _subgr-i of several subsets, N _total = ΣN _subgr-i Any one of them can select those showing “luciferase enzyme activity” that clearly shows a significant peak in a specific time zone. By determining similarity / commonality based on the presence of such characteristic peaks, the accuracy of “cluster analysis” described later can be increased.

  (5) Expression of “luciferase gene” of “reporter gene” based on temporal change of “luciferase enzyme activity” of each plant body determined to be “gene expression” in the “line” of interest "Cluster analysis" by pattern similarity

As a temporal change of “luciferase enzyme activity” of each plant body determined to be “gene expression” in the target “line”, time-series data (t _m [i + 1 ], P _{lum.corrected.-m} (t _m [i])), the temporal change of “luciferase enzyme activity” between each plant body (sample): S _Plant (m) in the line. Classification into multiple “clusters” based on similarity.

At that time, in each plant body (sample): S _Plant (m), the observed cell number N _{total of the} entire cell population is different, and in each individual cell, a “specific gene” and its allele When the expression of is progressed in parallel, the ratio of the expression frequency is not clear. In consideration of this point, time series data (t _m [i + 1], P _{lum.corrected.-m} (t _m [i]) that has been preliminarily selected and has been subjected to “measurement error correction” over a desired “long time” is selected. )
Measured light intensity: Numerical integration value of P _{lum.corrected.-m} (t _m [i])
∫P _{lum.corrected.-m} (t) dt
≈Σ1 / 2 · {P _{lum.corrected.-m} (t _m [i + 1]) + P _{lum.corrected.-m} (t _m [i])} · (t _m [i + 1] −t _m [i])
As a result, a time integral value is calculated. The value obtained by dividing the measured light intensity: P _{lum.corrected.-m} (t _m [i]) by this time integration value is an index indicating the temporal change of the normalized “luciferase enzyme activity”. . That is, during this “long time” measurement period, which time zone shows relatively high “luciferase enzyme activity”, that is, which time zone expresses the “specific gene”. Is an indicator of whether it is relatively high.

In which time zone, the number of cells expressing the “specific gene” is relatively high, that is, the similarity between the expression patterns of the “specific gene” in each plant (sample) Is based on time-series data indicating temporal changes in the normalized “luciferase enzyme activity”, {..., P _{lum.corrected.-m} (t _m [i]) / ∫P _{lum.corrected.-m} ( t) The “distance” between the vectors of dt,.

  Specifically, the calculation of the “distance” between vectors applies a row matrix calculation method.

As a result of separately testing the genotype of the plant grown from T1 seed in advance, it was confirmed that it showed a “heterotype” genotype having a chromosome configuration of (a ^* · a) or (a · a ^* ) With respect to the “line”, “cluster analysis” is performed on the plant that is grown from the T2 seed and is determined to have the expression of the “specific gene”.

At that time, as a sample of the "line", 2 ^four (16) a set, employing a total of two pairs, in each pair, by applying the techniques of various "cluster analysis", performed clustering after, together two sets, for a sample of 2 ⁵ (32), when subjected to the same clustered, in each set, the sample is judged to be most similar to the "pAIR", two sets of When combined, the ratio determined to be the most similar “pair” was also calculated.

  At that time, when the pair of samples that are determined to be the most similar “pair” in each group is combined, the ratio of the sample that is determined to be the most similar “pair” is also higher. The “analysis” method was selected as the optimal “cluster analysis” method.

As a method for calculating the “distance” between vectors,
The following four row matrix calculation methods:
・ Euclidian
・ Manhattan
・ Maximum
・ Canberra
Based on the calculated “distance” between vectors, as a coupling method used in the “clustering” process,
The following six bonding methods:
・ Average
・ Centroid
・ Complete
・ Mcquity
・ Single
・ Ward
Based on the above method, the optimal “cluster analysis” method combination was selected.

  As shown in Table 1 above, the combination using the Manhattan method as the inter-row matrix calculation method and the Ward method as the combination method is selected as the optimal “cluster analysis” method.

  Note that statistical analysis software R (http://www.r-project.org/) was used for “stratification” in cluster analysis based on the distance between each individual.

  In addition to the method selected above, the cluster analysis calculation method may include euclidean / maximum / cancella / binary / minkowski in the inter-row matrix calculation method and single / complete / average / mcquitity / in the combination method. median / centroid can be used as appropriate.

This “cluster analysis” method is shown in graph [B2] of FIG. 6. As a result of separately testing the genotype of the plant grown from T1 seed, (a ^* · a) or (a · a ^* ) Time of standardized “luciferase enzyme activity” of 16 plant bodies (samples) grown from the T2 seed for “line” that was confirmed to show a “heterotype” genotype having the chromosome structure of Applied to change data. As a result, as a result of “clustering”, hierarchical clustering shown in the dendrogram shown in the graph [B3] in FIG. 7 is performed. Finally, the thirteen plants (samples) determined to have “specific gene” expression are classified into two “groups”. The coupling distance between the two groups is 7.6 × 10 ⁵ as a threshold, and the plants (samples) classified into each group have substantially the same number of individuals.

The 13 plant bodies (samples) determined to have expression of this “specific gene” are the “homotype” of (a ^* · a ^* ) and (a ^* · a) or (a · a ^* ). It seems that it is classified into the genotype of “heterotype”.

As described above, for example, even when the genotype is a “heterotype” of (a ^* · a) or (a · a ^* ), a combination of alleles present on the homologous gene can be There is a “genotype”. Temporarily, regarding the “line”, all the plant bodies (samples) determined to have expression of the “specific gene” are “heterotype” of (a ^* · a) or (a · a ^* ). In some cases, these “genotypes” may cause “clustering” and eventually be classified into two “groups”. In that case, the coupling distance between the two groups is relatively short, and the threshold value is also relatively low.

In fact, in judging that it is significantly divided into two groups,
(I) Due to the hierarchical “clustering”, in the created dendrogram, the highest group is roughly divided into two groups with a threshold above a certain level.

  (Ii) At that time, the two groups roughly classified include at least two or more samples.

  It is necessary to satisfy the above two requirements.

For example, when one sample is classified into one group, statistically, a sample of “distance” below the threshold is set to 24 ⁴ (16) sets. The probability that it does not exist in the sample group is considerably low. On the other hand, if the upper threshold for dividing the top group does not reach a certain level, the probability that it can be judged that there is a significant difference between the two groups in statistics is considerably low.

(6) Within the target “line”, two types of “genetic expression” are roughly classified by “cluster analysis” based on temporal changes in “luciferase enzyme activity” of each plant body. Creation of “model waveform” representing the tendency of temporal change in “luciferase enzyme activity” in “Group” (6-1) “Unify time axis” operation About each plant body (sample) classified into each group The normalized time-dependent data of “luciferase enzyme activity” is (t _m [i], P _{lum.corrected.-m} (t _m [i]) / ∫P _{lum.corrected.-m} ( t) Time series data of dt), but for each plant (sample), the measurement time: t _m [i] is slightly different.

Therefore, for each plant body (sample), it is presumed that the normalized “luciferase enzyme activity” data was measured when it was measured at the same measurement time: t _LINE [i]. The data is converted into normalized time-dependent data of “luciferase enzyme activity” that “unify the time axis”.

Specifically, with respect to the unified time: t _LINE [i], the measurement time: t _m [i] in each plant body (sample) is around, and at least t _m [i -1] <t _LINE [i] <t _m [i + 1]. Focusing on this point, the standardized “luciferase enzyme activity” estimated at a unified time: t _LINE [i]: {P _{lum.corrected.-m} (t _LINE [i]) / ∫P _{lum.corrected.-m} (t) dt)} is estimated by interpolation based on the waveform data at times t _m [i−1], t _m [i], t _m [i + 1].

Specifically, the waveform data at time t _m [i−1], t _m [i], t _m [i + 1] monotonically increase in the waveform indicating the temporal change in the normalized “luciferase enzyme activity”. Alternatively, when monotonously decreasing or substantially constant, these three points are approximated by one straight line, and the value at time t _LINE [i] is employed as the estimated value on the approximate straight line. That is, by applying the least square method [method of least squares], the “approximate straight line” that linearly approximates three points on the waveform data at times t _m [i−1], t _m [i], and t _m [i + 1]. "

On the other hand, in the waveform showing the temporal change of the normalized “luciferase enzyme activity”, the waveform data at times t _m [i−1], t _m [i], and t _m [i + 1] are maximums of clear peaks. When it corresponds to a part, that is, when the time t _m [i] corresponds to the maximum point of the peak, it is estimated by interpolation using two points. Specifically, when t _m [i−1] <t _LINE [i] <t _m [i], a line connecting two points of time t _m [i−1] and t _m [i] The value at time t _LINE [i] is assumed to be an estimated value. When t _m [i] <t _LINE [i] <t _m [i + 1], the time t _LINE [i] is on a straight line connecting two points of time t _m [i] and t _m [i + 1]. The value at is the estimated value.

(6-2) “Model waveform creation” operation representative of the temporal change tendency of “luciferase enzyme activity” in two “groups” roughly classified by “cluster analysis”. From the data of time-dependent changes in “luciferase enzyme activity” that standardizes the “time axis” of each plant body (sample) group in each “line”, A “model waveform” representing the tendency of temporal change in “luciferase enzyme activity” is created.

Usually, the “waveforms” that are classified into the group and that have a similarity index “distance” that is substantially equal to the “plurality of waveforms” that are classified into the group are the “group”. “Model waveform” representing the commonality (trend) in “waveforms” belonging to That is, data obtained by chronologically arranging estimated normalized “luciferase enzyme activities” at a unified time: t _LINE [i]: {..., P _{lum.corrected.-m} (t _LINE [i ]) / ∫P _{lum.corrected.-m} (t) dt),...}}, The average of each group of plant bodies (samples) belonging to the group is generally determined by this “ This is a “model waveform” that represents the commonality (trend) among “waveforms” belonging to “group”.

However, the value of P _{lum.corrected.-m} (t _LINE [i]) is the original data: P _{lum.corrected.-m} (t _m [i]) itself has a statistical variation (variation). Therefore, it has the same level of statistical fluctuation (variation). For example, if the value of P _{lum.corrected.-m} (t _LINE [i]) has a Gaussian distribution (variation), the standard dispersion is {P _{lum.corrected.-m} (t _LINE [ i])} about ^1/2 . That is, the degree of relative _{variation, {P lum.corrected.-m (t} LINE [i])} 1/2 / P lum.corrected.-m (t LINE [i]) becomes, P _{Lum.Corrected The} smaller the value of _.-m (t _LINE [i]), the more it increases. Furthermore, when “normalization” is applied, the degree of relative variation in {P _{lum.corrected.-m} (t _LINE [i]) / ∫P _{lum.corrected.-m} (t) dt)} is As ∫P _{lum.corrected.-m} (t) dt) is smaller, the degree is further increased.

Considering this point, the simple “average” of the above-mentioned standardized “luciferase enzyme activity” is a standardized “plant” (sample) having a small ∫P _{lum.corrected.-m} (t) dt). It is easily affected by “variation” caused by “data”. Therefore, in the present invention, instead of the simple “average”, the “median” of “data” of each plant (sample) group belonging to the group is changed to the commonality (trend) in the “waveform” belonging to the “group”. ) As a representative value. This “median” is the normalized “luciferase enzyme activity” of each plant body (sample) at a unified time: t _LINE [i]: {P _{lum.corrected.-m} (t _LINE [i] ) / ∫P _{lum.corrected.-m} (t) dt)} are arranged according to their sizes and set as the central value. When the number of plants (samples) belonging to this group is 2n + 1, the (n + 1) th In the case of 2n, it is a simple “average” of the nth value and the (n + 1) th value.

For the group to which the number of plant bodies (samples) 7 in the dendrogram shown in the graph [B3] belongs to the graph [B4-1] shown in FIG. 8, the time unified according to the above procedure: t _LINE [i] In the above, “median value”: P _{lum.center-Gr1} (t _LINE [i]), time series representative data: (t _LINE [i], P _{lum.center-Gr1} (t _LINE [i] )) Is created. The graph [B4-1] illustrated in FIG. 8 also shows time-series data configured using a simple “average” for comparison.

Similarly, in the graph [B4-2] shown in FIG. 9, for the group to which the number of plant bodies (samples) 6 in the dendrogram shown in the graph [B2] belongs, the time unified according to the above procedure: t _LINE [ i], “median”: P _{lum.center-Gr1} (t _LINE [i]), time series representative data: (t _LINE [i], P _{lum.center-Gr1} (t _LINE [ i])) is created.

  The “model waveform” having the “median” as a representative value shows a maximum in the group in which the normalized “luciferase enzyme activity” of each plant body (sample) is commonly increased. It has become a thing. Actually, in this group, the time for which the normalized “luciferase enzyme activity” of each plant body (sample) takes a maximum is slightly around, and the maximum value is also dispersed. “Waveform” corresponds to “averaged” of these two variances.

  FIG. 10 shows a standardized “luciferase enzyme activity” in which all the plants (samples) included in the two groups included in the dendrogram shown in the graph [B2] are “unified time axis”. , And a time-series data composed of a “model waveform” and a simple “average” with “median” as a representative value for the first group shown in FIG. FIG. 9 also shows a “model waveform” having a “median value” as a representative value for the second group. At that time, the “model waveform” with “median” as a representative value for the first group and the “model waveform” with “median” as a representative value for the second group show the maximum values. Although the time zones are generally similar, it can be considered that there is a difference in the magnitude of the “maximum value” at each maximum peak.

  Focusing on the time zone showing the maximum, the regularity of the 24-hour period is estimated as a whole, and there are two peak pairs having a separation of about 12 hours in the 24-hour period. It is estimated to be. That is, each plant body (sample) to be measured is grown on the same 96-well plate and kept in the same environment, and the “regularity of a 24-hour period” similar to the so-called “circadian rhythm”. ", And the timing of the expression of the" specific gene "in the meantime.

  The created “model waveform” is stored in the model waveform table together with the result of “cluster analysis” indicating the configuration of the group, that is, the analysis result giving a tree diagram in “genotype automatic discrimination”.

In the process of creating the model waveform, in the “cluster analysis”, the temporal change in the “luciferase enzyme activity” of each plant body that is determined to have “gene expression” within the target “line”. Therefore, in each “line”, if the number of plants (samples) determined to have “gene expression” does not exceed at least 4, the two groups are roughly classified. Can not. In addition, if the criteria for determining that the two groups are roughly divided are not satisfied, the “model waveform” cannot be created for each group. In such a case, usually, for that "line", new, 2 ⁴ (16) was measured for a set of sample groups, based on the result, to analyze the same procedure.

  In the case where two independent sample groups are not divided into two groups in the same manner, the two sample groups are integrated and the same analysis proceeds. At that time, again, in the highest grouping, if one group belongs to only a single sample, this single sample is excluded, and “cluster analysis” is performed on the remaining sample group. It is determined whether or not a criterion for determining that the above two groups are roughly divided is satisfied. As a result, if it is determined that the remaining sample groups are roughly divided into the above two groups, a “model waveform” is created for these two groups.

  On the other hand, if it is not determined that the above two groups are roughly divided, it is determined that the “line” is a single group, and a “model waveform” is similarly created.

(7) In the target “line”, the “line” based on the “model waveform” representing the tendency of temporal change of “luciferase enzyme activity” in two “groups” roughly classified by “cluster analysis” Waveform analysis for the purpose of "feature grasp" (7-1) "Waveform decomposition" into a superposition of a plurality of "waveform functions" reflecting the waveform characteristics of "model waveform"
The “model waveform” created for each group includes information on a time zone in which “luciferase enzyme activity” commonly found in a plurality of plants (samples) classified into the group has a maximum. In addition, information on a representative “maximum value” in a time zone in which the commonly found “luciferase enzyme activity” exhibits a maximum is also included.

On the other hand, the temporal change in “luciferase enzyme activity” observed in the whole plant is originally the sum of the number of cells constituting the cell population N _total is the sum of the number of _sub- cells N _subgr-i , N _total = ΣN _subgr-i, and the cell group of the subset is a superposition of the temporal changes of the “luciferase enzyme activity” shown. Therefore, at first glance, the peak of the peak is crushed, and what can be regarded as a broad peak as a whole is inherently a multiple of temporal changes in the “luciferase enzyme activity” that shows the peak at different points in time. Can be interpreted.

  That is, the temporal change in “luciferase enzyme activity” in individual subsets of cells can be approximated by a unimodal waveform function exhibiting a similar half-width, and the unimodality exhibiting a similar half-width. As a result of a plurality of waveform functions overlapping, it can be interpreted that a waveform of a temporal change in “luciferase enzyme activity” observed in the whole plant is given.

In this individual subset of cells, the unimodal waveform function that approximately represents the temporal change in “luciferase enzyme activity” is a continuous function that is second-order differentiable, as illustrated above. For example, Lorentz function having peak height: h _peak , peak position: t _peak , full width at half maximum: 2Δt _half ; f (t) = h _peak × {1+ (t−t _peak ) ² / (Δt _half ) ² } ^{− 1}
And so on.

The created “model waveform” is in the form of time-series representative data: (t _LINE [i], P _{lum.center-Gr1} (t _LINE [i])). Although it reflects the tendency of change, it is not microscopically smoothed to the extent that numerical differentiation is possible. Originally, the temporal change of the standardized “luciferase enzyme activity” in individual plants is at least a continuous function that can be second-order differentiated, but the measurement error component at the time of measurement, that is, fine variation. Due to the included minute “noise”, the created “model waveform” also includes some minute “noise”.

  In consideration of these, at least when the numerical differentiation is performed, the waveform function has a continuity that can be second-order differentiated, so the moving average method is applied to the created "model waveform" to Apply smoothing. At that time, there are various methods in the moving average method that can suppress the distortion to the original waveform shape while removing the minute “noise” component from the waveform including the minute “noise” component. Although proposed, a polynomial fitting method is used here. At that time, at least when the numerical differentiation is performed, the number of smoothing points that can also achieve the objective of obtaining a waveform function having continuity that can be secondarily differentiated is examined. For example, the “model waveform” shown in FIGS. On the other hand, if 5 to 9 points are selected as the smoothing points, a minute “noise” component is removed, and at least when a numerical differentiation is performed, a waveform function having continuity capable of second order differentiation can be obtained. confirmed.

  Since the “model waveform” subjected to the smoothing process is smoothed by adopting a polynomial fitting method having a smoothing point of 5 or more, when it is numerically differentiated, it is in a state capable of third-order differentiation. . As described above, a smoothed “model waveform” can be approximated by superimposing a plurality of unimodal waveform functions that can be second-order differentiated, and a plurality of peaks of these unimodal waveform functions. It is estimated that a maximum peak or a shoulder peak is shown at a position corresponding to the position.

  At the maximum peak, the primary differential value changes from “positive → zero → negative”, and the secondary differential value shows a change of “large → small → large” type. In addition, even with a clear shoulder peak, the primary differential value changes as “positive → zero → negative → zero → positive” or “negative → zero → positive → zero → negative”. It shows a change from “small to large”. On the other hand, before and after a typical shoulder peak, the first-order differential value is a positive value indicating a change of “large → small → large” type, or a negative value indicating a change of “large → small → large” type. Therefore, the minimum point of the second derivative function and the time when the third derivative value becomes zero are specified as the time indicating the maximum peak or the shoulder peak of the “model waveform” subjected to the smoothing process.

An approximation of a “model waveform” that has been smoothed by applying the Davidson-Fletcher-Powell (DFP) method using the estimated time of the maximum peak or shoulder peak and the estimated value of the total number. As the waveform, the most suitable “composite waveform” is created. Here, the Lorentz function is adopted as the unimodal waveform function used to create this “synthetic waveform”, and each of the unimodal waveform functions has one unimodality with respect to the time showing the maximum peak or the shoulder peak. Corresponding to the maximum of the waveform function, the peak height: h _peak and the full width at half maximum: 2Δt _half of the unimodal waveform function at the peak position: t _peak are variously changed to create an optimal “composite waveform”. . At that time, optimization is performed so as to minimize the residual sum of squares between the “model waveform” subjected to the smoothing process and the “synthesized waveform”.

  The graph [B6-1] in FIG. 11 shows the result of smoothing the “model waveform” shown in FIG. 8 and the maximum peak in the “model waveform” after smoothing, or The result of creating an approximate waveform by “composite waveform” by superimposing a plurality of Lorentz functions in correspondence with the position of the shoulder peak is shown. The graph [B6-2] in FIG. 12 shows the result of smoothing the “model waveform” shown in FIG. 9 and the maximum peak in the “model waveform” after smoothing, or The result of creating an approximate waveform by “composite waveform” by superimposing a plurality of Lorentz functions in correspondence with the position of the shoulder peak is shown.

(7-2) Profile based on the result of “waveform decomposition” A plurality of “unimodal waveform functions”, that is, a “composite waveform” that is a superposition of Lorentz function type “decomposition curves” is smoothed. This corresponds to a “model waveform” that has been processed and “waveform decomposition”. That is, the characteristics of the “model waveform” subjected to the smoothing process, the total number of peaks, the position of each peak, the height, the half-value width, the integrated area of the entire “synthesized waveform”, and the smoothed process “ The integrated area of the “model waveform” is obtained as an analysis result.

  These analysis results and profiles are stored in the profile table 34. In addition, the graph display of this is stored in the file list table 32.

  The profile results shown in the graph [B6-1] in FIG. 11 and the graph [B6-2] in FIG. 12 show that the growth conditions, in particular, the light irradiation conditions, 12 days / 12 hours light / dark cycle 4 days ( 96 hours), and then in the measurement results under continuous light conditions for 3 days and a total of 7 days (168 hours), a 12-hour / 12-hour light-dark cycle of 4 days (96 hours) corresponds to the corresponding “circadian rhythm”. The peak of the 24-hour period is confirmed, and in addition, regularity of the 24-hour period is found even in the subsequent continuous light condition for 3 days.

(8) “Waveform comparison / classification” that verifies the similarity of temporal changes in the expression of the “specific gene” of interest between different “lines”
For each “line”, a different “line” is obtained by performing “cluster analysis” using the “model waveform” created from the measurement result of the temporal change in “luciferase enzyme activity” under each growth condition. It is also possible to verify the presence or absence of “similarity”.

  At that time, the “model waveform” of each line is divided by the integration area of the “model waveform” to obtain a standardized “model waveform”, and then “cluster analysis” in the above-mentioned “line”. Perform “cluster analysis” using the same procedure as above.

  An example in which “waveform comparison / classification” is performed on the graph [B7] shown in FIG. 13 to verify the similarity of temporal changes in the expression of the “specific gene” of interest between the different “lines”. Show. In order to show the degree of similarity, a standardized “model waveform” is shown in the form of “heat-map” along with its tree diagram.

  Hereinafter, the configuration and function of the analysis apparatus according to the present invention will be described in detail.

  In the specific example showing the embodiment of the present invention, not only the “waveform” corresponding to the temporal change in the “luciferase enzyme activity” is analyzed, but also “data” including the same measurement result and analysis result is obtained later. When a large number of data is accumulated, an index is added to the “data” for each “individual plant body” so that the analyzed information can be searched immediately. It is a system that can build the database that was done.

  Referring to FIG. 1, the present example is an input device 1 such as a keyboard, a data processing device 2 that operates under program control, a main storage device (database; DB) 3 that stores information, and an output of a display device or the like. The apparatus 4 is configured.

  The main storage device 3 includes a line table 31, a file list table 32, a waveform analysis list table 33, a profile table 34, and a model waveform table 35.

  The line table 31 issues an ID based on the line name of the line information 11 input from the input device 1. The table information includes information for specifying the measurement object of each measurement result, such as the line name, observation site, growth condition, phenotype, and the like.

  The file list table 32 collectively manages text files and image files created in each line. Based on the ID issued at the stage of creating the line table 31, a unique ID is held in the file list.

  The waveform analysis list table 33 combines the line table 31, the profile table 34, and the model waveform table 35 with a unique ID. In addition, it includes total area information and genotype information about the original waveform before waveform decomposition 20.

  The profile table 34 stores position information (time), height information (measured value), and half value width of each decomposition curve after the waveform decomposition 20. These pieces of information are managed based on the ID issued in the waveform analysis list table 33.

  The model waveform table 35 stores model waveform information (time, measured value) of each line after “model waveform creation” 19. These pieces of information are managed based on the ID issued in the waveform analysis list table 33.

  The computer (central processing unit; processor; data processing unit) 2 includes DB registration / update 12, data processing 6, expression level calculation 15, presence / absence of expression 16, model waveform creation 7, line characteristic grasping 8, graph creation 22 And waveform comparison / classification 23.

  Each of these means generally operates as follows.

  The DB registration / update 12 receives the data 10 output from the measuring device and the line information 11 separately input by the experimenter regarding the measurement object for each experiment. Next, the line table 31 in which the actual measurement “data 10” and “line information 11” in the experiment conducted in the past are stored is searched, and if there is existing information about the same “line”, Update. If there is no “line information 11” for “target line” as “line information 11”, new registration is performed.

  The stage of “data processing” 6 includes numerical data processing operations of “elapsed time conversion” 13 and “measurement error correction” 14, and “data 10” (time and measurement value = gene expression output from the measuring device). Process).

  Using the “graph creation” 22 function, a graph is created based on “primary processed data” obtained by performing the above-described processing on the actually measured “data 10” at the stage of “data processing” 6. .

  “Calculation of expression level” 15 determines the presence or absence of gene expression. The results (presence / absence of gene expression = phenotype) obtained for each plant (sample) to be measured are stored in the line table 31.

  The “presence / absence of expression” 16 branches depending on the presence / absence of gene expression determined in “calculation of expression level” 15. Only those recognized as “present” are subjected to the following “model waveform creation” 7 and “line feature grasp” 8.

  The stage of “model waveform creation” 7 includes numerical data processing operations of “genotype automatic discrimination” 17, “time axis unification” 18, and “model waveform creation” 19 based on plural types of waveform discrimination. First, in “genotype automatic discrimination” 17, it is judged by a plural type waveform discrimination operation whether or not a plurality of “expressed” samples (plants) included in a line show different waveforms depending on the genotype. To do. If it is determined that “shows a different waveform depending on the genotype”, a sample waveform is created for each type of waveform after classifying the sample (plant) for each “genotype”. If it is determined that “they do not show different waveforms depending on genotype” (they are the same), one model waveform is created. Here, the “model waveform” is a waveform that most reflects the characteristics of the line when a plurality of samples (plants) of the line draw different waveforms.

  The stage of “line feature grasping” 8 includes analytical operations of “waveform decomposition” 20 and “profile creation” 21, and features of the waveform based on the “model waveform” created in “model waveform creation” 7. Extract and create a “waveform profile”.

  In the stage of “waveform comparison / classification” 23, the model waveforms of the respective lines are compared and classified based on the data created in “line characteristic grasping” 8, in particular, “waveform profile”.

  Specific operations such as data exchange and storage in this analyzer will be described.

  Data [A1], data [B1] (part of actual data), and line information 11 (line: L00001, part: Seedling), which are data 10 output from the measuring instrument, provided from the input device 1. generation: T2, Condition: 4LD3Wc, Genotype: mixture, etc.) are passed to the DB registration / update 12.

  The DB registration / update 12 searches the DB based on this line information 11, updates if there is the same data, performs new registration if there is the same data, and issues an ID.

  Next, in “elapsed time conversion” 13, the time of the data output from the measuring device is formatted. The time display unique to the measuring device is unified into the HH: MM format (0 to 59 display is converted to 0 to 99 display for MM), and the measurement time is converted into the elapsed time from the measurement start time.

  In “measurement error correction” 14, automatic correction is performed when the measurement value output from the measuring instrument is incorrect. An instrument may detect an incorrect value due to various factors. Detect such false positive values (here, values that are considered to be expressed even though the gene is not expressed. Such values may occur due to static electricity generated during measurement), and correct the values. Do. The upper limit of the value to be corrected: i should be selected according to the target, and can be changed by “user designation” using the analysis system.

Conditions for correction:
(Formula 1-1) | X [n−1] −X [n + 1] | <20
(Formula 1-2) | (X [n] −X [n−1]) + (X [n] −X [n + 1]) | / 2> | i |
Value after correction:
(Formula 1-3) y = (X [n−1] + X [n + 1]) / 2

  As described above, before performing data processing 6 (data [A1], data [B1]) and after (data [A2], data [B2]), and after completing “measurement error correction” 14, sample data (data [A3]) is shown in graph [A3]. Upper limit of correction value: When i = 100, the data [A] is corrected to satisfy the conditional expression of “measurement error correction” 14. Data [B] does not satisfy the condition. For the purpose of explanation, the graph corresponding to the data before the measurement error correction 14 (data [A2], data [B2]) is converted into the graph [A2] and the graph [B2] and the data after the data (data [A3 ] A graph corresponding to [] is shown in graph [A3]. (Only the graph [B3], the graph [B5], the graph [B6], and the graph [B7] created in the waveform comparison / classification 23 created in the graph creation 22 are actually created in this system.) X axis of the graph Indicates time, Y-axis indicates gene activity, and each waveform indicates samples (1 to 16) included in one line. In order to make it easy to identify genes such as circadian rhythms, they are separated by broken lines and solid lines every 6 hours and 24 hours. Note that * in the data [A3] indicates that the outlier is corrected.

  By performing “measurement error correction” 14, an erroneous activity value (= false positive) can be eliminated. Further, it is possible to detect a waveform that has been overlooked as a relatively small value.

  Hereinafter, only the data [B] is handled in the embodiment of the operation after correction.

  The expression level calculation 15 determines the presence or absence of gene expression based on the measurement value and timing of the data [B3]. The condition for the presence or absence of expression varies depending on the user. For example, if the reference value “gene expression is present” is set to “when the measured value is 30 or more and 3 or more points are continuously observed”, the data [B3 ] Is “with gene expression” to satisfy this condition. The obtained result (presence / absence of gene expression = phenotype) is stored in the “line table” 31. By setting a threshold for the presence or absence of gene expression in the program, it is possible to prevent oversight or erroneous determination due to individual subjectivity.

  The “presence / absence of expression” 16 branches depending on the presence / absence of gene expression determined in “calculation of expression level” 15. Expression is considered to be present if at least one of the samples included in the line is recognized as “expressed”.

  In “multiple types of waveform determination” 17, it is determined whether or not a plurality of samples included in one line analysis draw different waveforms (gene activities) depending on the genotype. There are three types of genotypes: homo (2n): hetero (n): none, and when considered simply, gene expression (= gene activity, that is, a value that appears as a waveform) is half that of homo when the genotype is hetero. It becomes. However, in genes that self-suppress, gene expression is not simply halved, but a unique waveform may be drawn. Accordingly, first, the samples are classified into groups having the same genotype by using the function of “multiple types of waveform discrimination” 17, and “model waveform creation” 19 described below is executed for each sample.

Specifically, cluster analysis (*) was performed on samples (homo and hetero) in which activity was observed among the obtained samples,
1. When the dendrogram is cut at the specified threshold, it is divided into two groups.
2. Each group contains two or more samples,
When the above condition is satisfied, it is determined that there are different types of waveforms depending on the genotype. The condition of 1 is that if the groups are not separated when they are cut at a certain threshold (= distance, a measure of difference between groups), the degree of similarity between the groups cannot be said to be high to some extent. Similarly, when the group is divided into three or more groups, the degree of similarity between the groups is high, or there is a possibility of erroneous measurement or contamination. Condition 2 determines that there is a high probability of contamination if the number of samples obtained is one. Whether or not the waveforms differ depending on the genotype in the multiple-type waveform discrimination 17 is because the “shape” of the wave, that is, the timing of expression is more important than the expression level (those that greatly activate other genes in small amounts) In some cases, even if it is a large amount, it does not affect the living body very much), and the sample data is weighted and the total waveform area is all equal. This made it possible to put similar time-series gene expression in the same cluster with emphasis on timing. For cluster analysis, the calculation method varies depending on the purpose and application of the analysis, but as a result of applying cluster analysis to sample data whose genotype has already been determined by biological experiments (hereinafter, experimental data) Manhattan and the combination method are used because ward has the highest classification sensitivity and can accurately classify the genotype. The threshold value was determined empirically by subjecting the experimental data to cluster analysis and reasonably classifying it. The obtained result (the presence / absence of plural types of waveforms = phenotype) is stored in the “line table” 31. In this specific example, since the coupling distance between groups is 7.6e + 05, the group is divided into two groups, and each group includes two or more samples, so it is determined that the waveform differs depending on the genotype.

(*) Cluster analysis [cluster analysis]
Cluster analysis is a collective term for a method of collecting communities that are similar to each other among different objects that are similar to each other to create a community (cluster) and classifying the objects, and this is performed statistically. The result is expressed as a dendrogram (dendrogram) that combines similar data in sequence. Various methods (row matrix calculation method and combination method) have been proposed depending on the purpose and application of the analysis. The calculation method was determined and an analysis method using the statistical analysis software R was established. Cluster analysis is currently applied in various fields, such as classification of diseases based on symptoms and laboratory values, classification of companies by financial indicators, and classification of bacteria by shape and nature.

  In the “unify time axis” 18, the measurement time interval is unified between lines in the time formatted by the “elapsed time conversion” 13. In this case, similarities are compared between different lines later (waveform comparison / classification 23), so that different observation times are aligned between the lines in a pseudo manner. The least square method (*) was used to unify the time axis, minimizing errors before and after unification. The time width can be changed by “user specified”.

(*) Least squares method [method of least squares]
One of the approximate straight lines. A mathematical method of fitting a straight line that best represents the trend to a set of points (x, y). When n data (x ₁ , y ₁ ), (x ₂ , y ₂ ),..., (x _n , y _n ) are obtained, a straight line y = ax + b that best represents the tendency of the group is obtained. The given coefficients a and b can be expressed by Expression 2-1 and Expression 2-2.

  In the “model waveform creation” 19, a model waveform is created by using “unified time axis” data for each “genotype group” using the functions of “genotype discrimination” 18 and “unify time axis” 18. . From the measured values at the respective times after unifying the time axis of each sample, outliers (those that are abnormally separated from each other) that are recognized as “present” in the calculation of the expression level 15. If the same thing is analyzed, a certain value should be centered, but a distant value may be obtained.)

  As described above, the graph of the cluster analysis created in the “plurality waveform discrimination” 17 is the graph [B3], and the graph obtained by performing the “model waveform creation” 7 and adding the model waveform is the graphs [B4-1] and [B4-2]. Show. The graph [B3] shows a dendrogram (tree view) in which the horizontal axis indicates time, and the vertical axis indicates each sample and the similarity of the samples. In this cluster analysis, not only the dendrogram that connects each sample with similarity, but also when the gene was “up-regulated (activated)” or “down-regulated (inactivated)” at a glance. As can be seen, a plotting method that fuses pseudo images is used (heatmap). This is a relative color representation of the strength of each waveform, which indicates that the activity was stronger as blue → yellow → red. In this graph [B4], by creating a “model waveform”, it is possible to notice the “waveform properties” that are often overlooked. Further, in this example, it was found that the waveform engraved “circadian rhythm” every 12 hours. These waveform characteristics can be observed more prominently by “waveform decomposition” 20 described later.

  In the “waveform decomposition” 20, the waveform is decomposed and analyzed based on the “model waveform” created in the “model waveform creation” 19. When obtaining the characteristics (number of polar regions / position / height / half-value width / area, etc.) from the waveform, the waveform of such experimental data is not regular like a sound wave but has a complicated shape. It is difficult to take values with For example, even when taking the half-value width, the boundary of how far the width is from the polar region is ambiguous. Therefore, in order to express a complex waveform as a simpler waveform, it was considered to be a synthesis of a plurality of known curves and attempted to be adapted. In other words, one “model waveform” (complex curve) is decomposed into known functions (simple curves such as normal distribution) for each polar region.

  The “model waveform” created in “model waveform creation” 19 is used as the waveform to be decomposed. In the flow of the analysis method, first, the waveform is smoothed using the moving average method (this suppresses detection of erroneous polar regions). Next, the smoothed waveform is differentiated to identify the polar region. Based on the polar data, a known curve is adapted to each polar region using the DFP method which is one of the gradient methods. Finally, features (values) are acquired from these decomposed individual curves.

  Details of these flows will be described. First, since an uneven waveform corresponding to noise appears in the actual observed waveform, the value corresponding to noise is removed (corrected) using the moving average method. Although large noise has been removed in the measurement error correction 14, more detailed noise (such as that caused by the measurement environment and the measurement machine) is targeted. The moving average method as a denoising method is most commonly used and a plurality of methods have been proposed. Among them, a polynomial fitting method is used. The smoothing score varies depending on the target data, but 5 to 9 points showed good results. Next, the position and number of polar regions are specified by differentiating the waveform. A waveform obtained by removing noise by the moving average method (hereinafter referred to as an original waveform) is subjected to second and third differentiation. Since the minimum point of the second derivative curve and the zero position of the third derivative curve are the polar regions of the original waveform, the position and number are determined. Finally, it was assumed that the waveform could be expressed by a specific curve function (here, the Lorentz waveform based on the properties of the obtained waveform) for each polar region, and was fitted to this curve. The Lorentz waveform equation is shown in Equation 3. Here, ν represents the wave number, h represents the height of the polar region, u represents the wave number of the polar region, and w represents the half width. Based on the polar information obtained earlier, the Lorentz waveform is applied to the original waveform to calculate a Lorentz waveform (hereinafter, decomposition curve) that minimizes the error between the two curves.

  Finally, necessary information (number of polar regions / position / height / half-value width / area, etc.) is acquired from these decomposed individual decomposition curves. By performing “waveform decomposition” 20, each value relating to the extreme value can be acquired under the same conditions, and the extreme value can be clearly represented, so that the characteristic of the waveform that has been overlooked in the past can be found.

  In the profile creation 21, the characteristics for each line extracted by the waveform decomposition 20 are stored in the profile table 34 based on the line table 31 as a profile. Currently, the data stored as a profile includes ID, the number of polar regions, the position of the polar region of the original waveform (time), the value of the polar region of the original waveform, the position of the polar region of the decomposition curve (time), and the polar region of the decomposition curve , Half width, area of original waveform, total area of decomposition curve (= composite waveform of decomposition curve, sum of decomposition curves).

  The graph creation 22 receives the above results and creates a graph. In addition to a graph of the data after the data processing 6 is completed, a model waveform and a decomposition curve are added for a plurality of waveforms when the plural types of waveform discrimination 17 are recognized, and for the original waveform otherwise. In the case of a line in which expression is not observed due to presence / absence of expression 16, only a graph of data after data processing 6 is created. The display method was divided at regular intervals and processed so that it was easy to confirm growth conditions and observe gene expression with a circadian rhythm. The created graph is stored in the file list table 32 based on the already registered line table 31 information.

  As described above, the graphs after “waveform decomposition” 20 are shown in graphs [B6-1] and [6-2] (there are two types of “multiple types of waveform discrimination” 17 because the waveforms differ depending on “genotype”). X axis is time, Y axis is "gene activity", Original (red or blue thick line) of the waveform is the model waveform, Separate ~ is the decomposition curve of the model waveform, Sum is the synthesis curve again (Sum of each decomposition curve) is shown. The accuracy of the decomposition curve is better as the waveforms of Original and Separate are more similar. By performing “waveform decomposition” 20, not only can the half width be accurately obtained, but in this example, it can be seen that “circadian rhythm” is carved even after 96 hours, which was conventionally overlooked.

  “Waveform comparison / classification” 23 performs cluster analysis between different lines, and compares and classifies similarities. The lines to be compared can be selected according to conditions using the information in the line table 31 (“light conditions apply to light / dark cycle 4 days + continuous light conditions 3 days”, etc.). Specifically, using the line table 31 and the model waveform table 35, “model data” that matches each line condition is extracted. Similar to “Multiple Waveform Discrimination” 17, when comparing between lines, the timing of expression is more important than the expression level. Therefore, the “model data” is weighted so that the total area of the waveforms is all equal. Data "was the subject of analysis. The “cluster analysis” is performed in the same manner as “multiple types of waveform discrimination” 17 using a heat map and the calculation method of the cluster analysis.

  The graph after the “waveform comparison / classification” 23 is shown in the graph [B7]. The horizontal axis represents time, and the vertical axis represents “model data” of each line and a dendrogram (dendrogram) showing the similarity of samples. Thereby, the relationship between each gene can be estimated more easily. In this example, the gene can be broadly divided into a gene “activated once for a certain period and then inactivated” and a gene “repeatedly activated and deactivated at certain intervals”. If you look more closely, it may be possible to find a gene that promotes activation / inactivation of a gene from the time difference.

  The “waveform comparison / classification” 24 branches depending on whether a normal search is performed or the process proceeds to the waveform comparison / classification 23. When performing a normal search, the result processed in the waveform analysis 5 described above is searched from the DB. The search target is all of the line table 31, the file list table 32, the waveform analysis list table 33, the profile table 34, and the model waveform table 35. Information provided from the outline to the details depending on "user needs" Can be selected. It is also possible to narrow down the lines to be browsed according to conditions. When performing waveform comparison / classification, the process proceeds to waveform comparison / classification 23.

  “Search” 25 selects contents to be processed. Select whether to execute “search” or “waveform comparison / classification” 23 and set the detailed conditions.

  “Screen output” 26 presents information analyzed so far. The presentation line, details of information, acquisition of graphs, etc. can all be selected by the user.

  Below, explanation is added about the advantage by utilizing the analysis system concerning the present invention.

In the conventional analysis method, the numerical data obtained from the experimental measuring instrument manually can be arranged and calculated by a spreadsheet software such as Excel, and only a graph can be created. When features are extracted manually from the graph, visual errors and omissions occur, and the subjectivity of each worker is included in the data. Moreover, these data are obtained in large quantities, and with these data,
・ It takes time ・ It takes labor ・ An error occurs in the numerical value ・ Subjectivity of the measurer enters ・ There is a lack of detection of graph features
・ Understanding the characteristics of each line is difficult (or impossible)
・ Comparison between lines is difficult (or impossible)
・ It is difficult (or impossible) to process a large amount of experimental data
There was a problem.

In the analysis system according to the present invention, by systematizing all the above functions, for the problem in the conventional analysis method,
・ It takes time → Significantly reduces ・ It takes labor → Significantly decreases ・ Errors in numerical values → Significantly decreases ・ Increased subjectivity of the measurer → All objectives ・ Detection omission → No detection omission → All numerical (quantitative)
And the problem can be solved collectively.

  In addition, the characteristics of each line can be quantitatively grasped from the profile, thereby enabling comparison and classification between the lines, and inferring the relationship between genes that cause these gene expression patterns.

Moreover, the following advantages can be obtained by applying the analysis system according to the present invention.
・ When a large amount of time-series expression data is obtained, the lines can be compared and classified with quantitative data that has been objectively analyzed.
・ By performing waveform decomposition, a single waveform is viewed as a plurality of waveforms, so that the characteristics can be seen more remarkably. In fact, results such as the discovery of “circadian gene expression” that could not be found by conventional methods were observed.
-It can also be used when observing gene expression over time in a microarray or cell array.
・ Estimates the order and mode of expression between genes and contributes to life science.

  The analysis method of the present invention can be applied to detailed statistical analysis for the purpose of verifying the similarity of a target for drawing a continuous waveform, such as time-series gene expression analysis, regardless of the experimental method.

It is a figure which shows typically the structure of the analysis system concerning this invention. " In the analysis method of the present invention, the “data 10” (time and time) output from the measuring instrument by the numerical data processing operation of “elapsed time conversion” 13 and “measurement error correction” 14 in the stage of “data processing” 6. FIG. 4 is a diagram for explaining the processing of A1 → A2 → A3. In the analysis method of the present invention, “data 10” (B1) processed by numerical data processing operations of “elapsed time conversion” 13 and “measurement error correction” 14 at the stage of “data processing” 6 and post-processing FIG. 4 is a diagram for explaining “primary processing data” (B2) used for “calculation of expression level” 15; “Data 10” output from the measuring device is a graph [A2] in which the time change of “luciferase enzyme activity” of each sample after the processing of “elapsed time conversion” 13 is also displayed in a graph. The “primary processing” of each sample processed from the “data 10” output from the measuring instrument by the numerical data processing operation of “elapsed time conversion” 13 and “measurement error correction” 14 in the stage of “data processing” 6. Graph [A3] in which “data” is also displayed as a graph is shown. A graph [B2] is also shown, in which a temporal change of the “luciferase enzyme activity” normalized by the integral value of the “luciferase enzyme activity” of each sample used for “cluster analysis” is also displayed as a graph. The graph [B3], in which the dendrogram created by “cluster analysis” according to the classification and the temporal change in the normalized “luciferase enzyme activity” in each sample, is also displayed in a heatmap format. Graph [B4-1], in which the “model waveform” created for the sample group belonging to the upper group among the samples roughly divided into two groups by “cluster analysis” shown in graph [B3] is displayed as a graph. Indicates. Graph [B4-2], in which the “model waveform” created for the sample group belonging to the lower group among the samples roughly divided into two groups by “cluster analysis” shown in graph [B3] is displayed as a graph. Indicates. The “model waveform” created for the sample group belonging to the upper group showing the graph [B4-1] and the “model” created for the sample group belonging to the lower group showing the graph [B4-2] The graph [B5] is displayed in a graph by comparing with the “waveform”. For the sample group belonging to the upper group, which represents the graph [B4-1], a “model waveform” that has been smoothed based on the created “model waveform”, and a plurality of “decomposition curves” that have been “waveform decomposed” ”, A graph [B6-1] in which the“ composite waveform ”obtained by synthesizing the“ decomposition curve ”is compared and displayed as a graph. For the sample group belonging to the lower group, which represents the graph [B4-2], a “model waveform” that has been smoothed based on the created “model waveform”, and a plurality of “decomposition curves” that have been subjected to “waveform decomposition” , And a graph [B6-2] in which the “composite waveform” obtained by synthesizing the “decomposition curve” is displayed as a graph. Using the “model waveform” created for each “line” and showing the temporal change in “luciferase enzyme activity”, “waveform comparison / classification” is performed between multiple “lines”. A graph [B7] in which the dendrogram showing the result and the “model waveform” are displayed together in a heatmap format is shown.

Claims (8)

An analysis method for the purpose of comparative classification based on feature extraction of time series gene expression data in individuals derived from the same "line" and characteristics of time series gene expression data between individuals,
The time-series gene expression level data in the individual organism uses a luciferase gene as a `` reporter gene '' capable of non-destructively monitoring the expression level of the gene, and as a chronological change in `` luciferase enzyme activity '', Observed data,
The analysis method is
It is determined whether or not “measurement error” data caused by the measuring instrument is mixed in the time-series numerical data of “luciferase enzyme activity” obtained from the measuring instrument. A process of “measurement error correction”, which performs an operation of correcting the “measurement error” data to an “estimated value” based on the numerical data before and after that;
Comparing numerical data of time-series changes of “luciferase enzyme activity” for multiple organisms originating from the same “line”, hierarchical clustering is performed according to the level of similarity, and at least Performing “cluster analysis” operation to form one or two groups including the individual of “multiple types of waveform discrimination”;
Based on numerical data of time-series changes in the “luciferase enzyme activity” for biological individuals originating from the same “line”, observed in the individual at a desired elapsed time starting from the measurement start time It is converted into numerical data of time-series changes in predicted “luciferase enzyme activity”, and “time axis unification” is performed between the numerical data of time-series changes in “luciferase enzyme activity” in each individual. Process of “unification of time axis”;
Using the numerical data of the time-series changes in “luciferase enzyme activity” in each individual for which “time axis unification” has been made, a plurality of types determined to belong to the same group by the “plurality waveform discrimination” Create a “model waveform” that shows the time-series changes in the representative value of “luciferase enzyme activity”, reflecting the similarity of the time-series changes in “luciferase enzyme activity” in this group. The process of “model waveform creation”;
Based on a “model waveform” indicating a time-series change of a representative value of “luciferase enzyme activity” for a group consisting of a plurality of individuals determined to belong to the same group, a plurality of “unimodal waveform functions” are provided. Overlay, create a “composite waveform” that approximately shows the waveform characteristics of the “model waveform”, and for each of the “unimodal waveform functions”, the values of the peak position, peak height, and half width A method of analyzing time-series gene expression level data, comprising a step of “waveform decomposition” to determine Using the “model waveform” created in the process of “model waveform creation” for the biological group derived from each “line”,
Compare “model waveforms” of organisms derived from multiple “lines”, perform hierarchical clustering according to the level of similarity, and at least one group that includes multiple “lines”, or The method of analyzing time-series gene expression level data according to claim 1, further comprising a “waveform comparison / classification” step of performing two “cluster analysis” operations. The analysis of time-series gene expression level data according to claim 1, wherein the "unimodal waveform function" used in the "waveform decomposition" step is a Lorentz function type waveform function. Method. In the process of “model waveform creation”,
In a group consisting of a plurality of individuals discriminated to belong to the same group, as the representative value of “luciferase enzyme activity” that reflects the similarity of time-series changes in “luciferase enzyme activity” in this individual group, the individual 2. The method for analyzing time-series gene expression level data according to claim 1, wherein a central value of “luciferase enzyme activity” at each time in a group is selected. In the “Multiple Waveform Discrimination” process,
When comparing numerical data of time-series changes in the “luciferase enzyme activity” for multiple organisms that originate from the same “line”, each individual can be clustered according to the level of similarity. The similarity of numerical data of time-series changes in “luciferase enzyme activity” between the two is calculated using the inter-row matrix calculation method to calculate the “distance” between the numerical data of time-series changes in “luciferase enzyme activity”. The time-series gene expression level data analysis method according to claim 1, wherein the time-series gene expression level data is used. As the inter-row matrix calculation method, a manhattan method is used,
The time series according to claim 6, characterized in that the Ward method is used as a binding method for clustering based on the "distance" between the numerical data of the time series changes of the "luciferase enzyme activity". Analysis method of gene expression level data. An analysis system that can be used for analysis for the purpose of comparative classification based on feature extraction of time-series gene expression data in organism individuals originating from the same "line" and characteristics of time-series gene expression data between individuals Because
The analysis system, as a mechanism for performing analysis according to the time series gene expression level data analysis method according to claim 1,
It is determined whether or not “measurement error” data caused by the measuring instrument is mixed in the time-series numerical data of “luciferase enzyme activity” obtained from the measuring instrument. A “measurement error correction” mechanism that performs an operation to correct the “measurement error” data to an “estimated value” based on the numerical data before and after that;
Comparing numerical data of time-series changes of “luciferase enzyme activity” for multiple organisms originating from the same “line”, hierarchical clustering is performed according to the level of similarity, and at least A "multiple-type waveform discrimination" mechanism that performs a "cluster analysis" operation to form one or two groups, including
Based on numerical data of time-series changes in the “luciferase enzyme activity” for biological individuals originating from the same “line”, observed in the individual at a desired elapsed time starting from the measurement start time It is converted into numerical data of time-series changes in predicted “luciferase enzyme activity”, and “time axis unification” is performed between the numerical data of time-series changes in “luciferase enzyme activity” in each individual. "Unified time axis"mechanism;
Using the numerical data of the time-series changes in “luciferase enzyme activity” in each individual for which “time axis unification” has been made, a plurality of types determined to belong to the same group by the “plurality waveform discrimination” Create a “model waveform” that shows the time-series changes in the representative value of “luciferase enzyme activity”, reflecting the similarity of the time-series changes in “luciferase enzyme activity” in this group. The “model waveform creation” mechanism;
Based on a “model waveform” indicating a time-series change of a representative value of “luciferase enzyme activity” for a group consisting of a plurality of individuals determined to belong to the same group, a plurality of “unimodal waveform functions” Overlay, create a “synthetic waveform” that approximately shows the waveform characteristics of the “model waveform”, and for each of the “unimodal waveform function” values for the peak position, peak height, and half width A system for analyzing time-series gene expression data, characterized by having a “waveform decomposition” mechanism for determining Using the “model waveform” created in the process of “model waveform creation” for the biological group derived from each “line”,
Compare “model waveforms” of organisms derived from multiple “lines”, perform hierarchical clustering according to the level of similarity, and at least one group that includes multiple “lines”, or The time series gene expression level analysis system according to claim 7, further comprising a “waveform comparison / classification” mechanism for performing two “cluster analysis” operations.