JP2010108402A - Device for analysis of two types of time-series data, and computer-readable storage medium recording program for analysis of two types of time-series data - Google Patents
Device for analysis of two types of time-series data, and computer-readable storage medium recording program for analysis of two types of time-series data Download PDFInfo
- Publication number
- JP2010108402A JP2010108402A JP2008282192A JP2008282192A JP2010108402A JP 2010108402 A JP2010108402 A JP 2010108402A JP 2008282192 A JP2008282192 A JP 2008282192A JP 2008282192 A JP2008282192 A JP 2008282192A JP 2010108402 A JP2010108402 A JP 2010108402A
- Authority
- JP
- Japan
- Prior art keywords
- analysis
- data
- segment
- series data
- time series
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
Images
Landscapes
- Complex Calculations (AREA)
Abstract
Description
本発明は、2系列の時系列データの解析装置及び2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体に関するものである。 The present invention relates to a two-series time-series data analysis apparatus and a computer-readable recording medium that records a two-series time-series data analysis program.
本発明者は、先に、時系列データを高精度に解析可能な解析装置として、例えば特許文献1、特許文献2に示されるように、最大エントロピー法(MEM)を利用した時系列データの解析手法及びその解析手法を適用したコンピュータプログラムを実行する解析装置を提案している。尚、最大エントロピー法を利用した時系列データの解析手法は、例えば非特許文献1、非特許文献2にも記載されている。
The inventor previously analyzed time-series data using a maximum entropy method (MEM) as shown in
上述した解析手法は、解析対象の時系列データを一般化三角多項式により表現して解析を行うもので、概ね、時系列データに対して最大エントロピー法を適用して、高精度なパワースペクトル密度(PSD)を求め、その顕著なピークから、一般化三角多項式を解く初期値としての、項数と、各三角項の周期値を求める過程と、これらの初期値により、線形化された非線形最小自乗法を用いて、各三角項のパラメータを求めると共に、残差時系列データを求め、これらを加算して時系列データを再構成する過程とから構成されている。 The analysis method described above is to analyze the time series data to be analyzed by expressing it with a generalized trigonometric polynomial. Generally, the maximum entropy method is applied to the time series data, and the power spectral density ( PSD), and the number of terms as the initial value for solving the generalized trigonometric polynomial and the periodic value of each triangular term from the prominent peak, and the nonlinear minimum self-linearized by these initial values. The method includes a process of obtaining parameters of each triangular term by using multiplication, obtaining residual time series data, and adding these to reconstruct time series data.
この解析手法及びそれを適用した解析装置では、何らの数値モデルも想定せず、時系列データのみを用いて高精度の解析を行え、しかもその解析は、通常のパーソナルコンピュータ程度の演算能力でも十分に実用的であるという特徴を有している。
しかしながら以上の文献に記載されている解析装置は、一系列の時系列データを解析するのに用いられるものであり、本発明者は、これらの解析装置に用いられている解析手法の上述した特徴を更に発揮させるべく、時系列データのみを用い、時系列に関わるその他の情報を一切使用することなく、二つの時系列の相互相関にかかわる知見を得ることのできる解析装置を提供することを課題とした。 However, the analysis devices described in the above documents are used for analyzing a series of time-series data, and the inventor has described the above-described features of the analysis methods used in these analysis devices. To provide an analysis device that can obtain knowledge related to the cross-correlation of two time series without using any other information related to the time series. It was.
ところで上述した解析装置で一系列の時系列データの解析を行う場合、演算速度、演算精度、記憶容量等の点でコンピュータの能力に限界があるため、長大な時系列データを一括して処理するのが困難であり、この場合には、長大な時系列データを分割して処理する必要がある。一方、時系列データには、観測期間内のデータ全体に渡ってデータ構造が均一に保たれるものの他、観測期間内でデータ構造が時々刻々と変化するものがあり、後者の場合には、観測期間内のデータ全体を一括して処理するのは困難である。 By the way, when analyzing a series of time series data with the above-described analysis device, the computer's ability is limited in terms of computation speed, computation accuracy, storage capacity, etc., so long time series data is processed in a lump. In this case, it is necessary to divide and process long time-series data. On the other hand, time-series data includes data whose data structure is kept uniform over the entire data in the observation period, as well as those whose data structure changes from moment to moment within the observation period. In the latter case, It is difficult to process all the data in the observation period at once.
そこで本発明では、解析対象の時系列データを、そのデータ点数やデータ構造に応じて適切にセグメント化して処理を行うことにより、解析する時系列データによらず、常に適切に高精度な解析を行えるようにすることも課題とした。 Therefore, in the present invention, time series data to be analyzed is appropriately segmented according to the number of data points and the data structure and processed, so that it is always possible to appropriately and accurately analyze regardless of the time series data to be analyzed. Making it possible was also an issue.
尚、本発明が対象とする時系列データとは、実在の時系列データ、即ち、ある限られた期間に観測されたある点数の有限長離散時系列のデータであり、理論上の存在である無限長のデータや連続データは対象としない。また本発明が対象とする時系列データは、上述したような通常の時系列データの他、例えば地層断面の距離に対する輝度変化等の、空間的に連続するデータも対象とするものである。 The time-series data targeted by the present invention is real time-series data, that is, data of a finite-length discrete time series with a certain number of points observed in a certain limited period, and is theoretically present. Infinite length data and continuous data are not included. The time-series data targeted by the present invention is not limited to normal time-series data as described above, but also includes data that is spatially continuous, such as a change in luminance with respect to the distance of the formation cross section.
上述した課題を解決するために、本発明では、解析対象の2系列の時系列データを、夫々解析して、その一般化三角多項式パラメータを求める2系統の単一系列解析手段と、夫々の単一系列解析手段により求められた2系列の時系列データの一般化三角多項式パラメータから、2系列の時系列データの間の相互相関を定量化する諸量を求める相互相関解析手段とから構成し、単一系列解析手段は、解析対象の時系列データを解析するシーケンスを決定する解析シーケンス決定手段と、解析シーケンスに従い、最大エントロピー法及び非線形最小自乗法によりデータを解析して一般化三角多項式のパラメータと残差を求めるデータ解析手段とから構成し、相互相関解析手段は、一般化三角多項式パラメータからパワースペクトル密度を求めるパワースペクトル密度計算手段と、クロススペクトル(コスペクトル、クオドスペクトル)とコヒーレンスを含む相互相関を定量化する諸量を求める相互相関諸量計算手段と、相互相関諸量計算手段により求められたコヒーレンスを理論値と比較して、それらの差が所定値以下の場合に、その差をなくすようにクロススペクトルを補正する整合性評価補正手段と、整合性評価補正手段において評価した差が所定値以上の場合に、単一系列解析手段に再解析の指令を出力し、所定値以下の場合には、相互相関諸量計算手段により計算された諸量と、パワースペクトル密度計算手段により計算されたパワースペクトル密度及び整合性評価補正手段により補正されたクロススペクトルを出力する分岐出力手段を構成した2系列の時系列データの解析装置を提案する。 In order to solve the above-described problems, the present invention analyzes two series of time series data to be analyzed, respectively, and obtains a generalized trigonometric polynomial parameter of two systems, and each of the single series analysis means. A cross-correlation analysis means for obtaining various quantities for quantifying the cross-correlation between the two series of time series data from the generalized trigonometric polynomial parameters of the two series of time series data obtained by the one series analysis means; The single series analysis means includes an analysis sequence determination means for determining a sequence for analyzing the time series data to be analyzed, and according to the analysis sequence, the data is analyzed by the maximum entropy method and the nonlinear least squares method, and the parameters of the generalized trigonometric polynomial And a data analysis means for obtaining the residual, and the cross correlation analysis means is a parameter for obtaining the power spectral density from the generalized trigonometric polynomial parameters. -Theory of spectral density calculation means, cross-correlation quantity calculation means for quantifying cross-correlation including cross spectrum (cospectrum, quad spectrum) and coherence, and coherence obtained by cross-correlation quantity calculation means When the difference between these values is less than or equal to a predetermined value, the consistency evaluation correction means for correcting the cross spectrum so as to eliminate the difference, and the difference evaluated by the consistency evaluation correction means is greater than or equal to the predetermined value If the re-analysis command is output to the single series analysis means and the value is less than the predetermined value, the quantities calculated by the cross-correlation quantity calculation means and the power spectrum density calculated by the power spectrum density calculation means And a two-series time-series data analysis device comprising branch output means for outputting a cross spectrum corrected by the consistency evaluation correction means Suggest.
また本発明では、上記の解析装置において、パワースペクトル密度計算手段は一般化三角多項式パラメータからパワースペクトル密度を求める計算において、一般化三角多項式の夫々の項の夫々につき、δ函数の原型を用いて全周波数帯に裾の広がる夫々一つのスペクトルピークを構成することにより、パラメータから得られる各周波数に対応する孤立ピーク群を構成し、その重畳として時系列データのスペクトルを求める構成とする2系列の時系列データの解析装置を提案する。 In the present invention, in the above analysis device, the power spectral density calculation means uses the prototype of the δ function for each term of the generalized trigonometric polynomial in the calculation for obtaining the power spectral density from the generalized trigonometric polynomial parameters. By configuring each spectrum peak that spreads in the entire frequency band, an isolated peak group corresponding to each frequency obtained from the parameters is formed, and two series of spectra are configured to obtain the spectrum of the time series data as its superposition. A time-series data analysis device is proposed.
また本発明では、上記の解析装置において、δ函数の原型として、余弦函数に微弱な雑音を加えた時系列をラグ2にて最大エントロピー法のスペクトル密度の計算を行って得られたスペクトルを利用する2系列の時系列データの解析装置を提案する。
Further, in the present invention, in the above analysis apparatus, a spectrum obtained by calculating the spectral density of the maximum entropy method at
また本発明では、上記の解析装置において、解析シーケンス決定手段は、解析対象の時系列データを、予め設定された最大データ点数を越えないようにセグメントに分割する分割処理過程と、各セグメントに対してデータ構造の均一性を評価する評価処理過程と、この評価に基づき、不均一のセグメントデータを均一の複数のセグメントに分割する分割処理過程と、分割された均一のセグメントの夫々を基本セグメントとして、それらに対して、記憶手段に記憶領域を設定する記憶領域設定処理過程を備えている2系列の時系列データの解析装置を提案する。 According to the present invention, in the above analysis device, the analysis sequence determination means divides the time series data to be analyzed into segments so as not to exceed a preset maximum number of data points, and for each segment, The evaluation process for evaluating the uniformity of the data structure, the division process for dividing the non-uniform segment data into a plurality of uniform segments based on this evaluation, and each of the divided uniform segments as a basic segment In response to this, a two-series time-series data analysis device having a storage area setting process for setting a storage area in the storage means is proposed.
また本発明では、上記の解析装置において、セグメントは、解析対象の時系列データの全期間データに対応する最上位レイヤーのセグメントから、順次2分割されて構成される、下位のレイヤーの夫々に含まれる2の累乗のセグメントと、最下位のレイヤーにおいて、データ構造が変化するセグメントを、データ構造の変化している個所において分割して順次形成されるサブレイヤーの各セグメントとから構成される2系列の時系列データの解析装置を提案する。 According to the present invention, in the above analysis apparatus, the segment is included in each of the lower layers configured by sequentially dividing the segment from the highest layer segment corresponding to the entire period data of the time series data to be analyzed. 2 series comprising segments of powers of 2 and sub-layer segments sequentially formed by dividing a segment whose data structure changes in the lowest layer at a location where the data structure changes We propose a time series data analysis device.
また本発明では、上記の解析装置において、データ構造の均一性は、下式による最大ラグ値までの最大エントロピー法のスペクトル密度の計算におけるバーグの係数Pmの、ラグ値に対する推移により評価する2系列の時系列データの解析装置を提案する。
また本発明では、上記の解析装置において、解析シーケンス決定手段には、不等間隔の時系列データを等間隔化処理する等間隔化処理過程を設けた2系列の時系列データの解析装置を提案する。 Further, in the present invention, in the above analysis apparatus, a two-series time-series data analysis apparatus is proposed in which the analysis sequence determination means is provided with an equal-interval processing process for equal-interval processing of non-uniformly-interval time-series data. To do.
また本発明では、上記の解析装置において、データ解析手段は、解析対象の時系列データから取り出したセグメントデータに対して下式による最大ラグ値までの最大エントロピー法のスペクトル密度の計算におけるバーグの係数Pmの、ラグ値に対する推移を求める計算過程と、求めた推移から最大エントロピー法によるパワースペクトル密度の計算に用いるラグ値を決定するラグ値決定過程と、決定されたラグ値により上記パワースペクトル密度を求めるパワースペクトル密度計算過程と、求められたパワースペクトル密度から主要ピークの周波数と、そのパワーを抽出する抽出過程と、抽出された周波数とパワーとから非線形最小自乗法により当てはめ曲線のパラメータを求めるパラメータ導出過程と、求められた当てはめ曲線とセグメントデータとを比較して一致度を評価し、所定の一致度以上の場合には上記当てはめ曲線のパラメータをセグメントに対応する記憶手段に記憶して次の過程に移行し、所定の一致度に達しない場合には上記当てはめ曲線のパラメータを破棄して上記ラグ値を決定する過程に移行する評価過程とを備えており、ラグ値決定過程は、設定される解析の精度に応じて、より高精度の解析において、より大きいラグ値を選択するように構成し、抽出過程は、設定される解析の精度に応じて、より高精度の解析において、より多くのピークを主要ピークとして抽出するように構成すると共に、評価過程は、設定される解析の精度に応じて、より高精度の解析において、より高い一致度となるように構成し、相互相関解析手段からの再解析指令により、高精度の解析が設定されるように構成した2系列の時系列データの解析装置を提案する。 Further, in the present invention, in the above analysis apparatus, the data analysis means includes the Burg coefficient in the calculation of the spectral density of the maximum entropy method up to the maximum lag value according to the following equation for the segment data extracted from the time series data to be analyzed: The calculation process for determining the transition of P m with respect to the lag value, the lag value determining process for determining the lag value used for calculating the power spectral density by the maximum entropy method from the determined transition, and the power spectral density based on the determined lag value Calculate the parameters of the fitting curve by the nonlinear least square method from the power spectral density calculation process to find the frequency of the main peak from the obtained power spectral density, the extraction process to extract the power, and the extracted frequency and power The parameter derivation process and the obtained fitting curve and segment The degree of coincidence is evaluated by comparing with the data, and if the degree of coincidence is equal to or higher than the predetermined coincidence, the parameters of the fitting curve are stored in the storage means corresponding to the segment, and the process proceeds to the next process, and the predetermined coincidence is reached If not, the evaluation curve transitions to the process of determining the lag value by discarding the parameters of the fitting curve, and the lag value determination process has a higher accuracy depending on the accuracy of the set analysis. The analysis process is configured to select a larger lag value, and the extraction process is configured to extract more peaks as the main peak in the higher-accuracy analysis according to the set analysis accuracy. In addition, the evaluation process is configured to have a higher degree of coincidence in the higher-accuracy analysis according to the set accuracy of analysis, and a high-precision analysis is performed by a re-analysis command from the cross-correlation analysis means. Analysis suggest analyzer of the time-series data of the configuration the two series to be set in.
また本発明では、上記の解析装置において、データ解析手段におけるセグメントデータの解析は、最上位のレイヤーから最下位のレイヤーに向かって順次行うものとし、この際、上位のレイヤーのセグメントのセグメントデータは、セグメント当たりのデータ点数が予め設定された最大データ点数を越えないように、解析対象の時系列データのデータ点数をバンチング処理により低減してセグメントデータを構成する処理を行う構成とすると共に、各レイヤーの各セグメントに対応して一時記憶領域を構成し、上位側のセグメントに対して求められた多項式パラメータを下位の各レイヤーの各セグメントの一時記憶領域に記憶すると共に、上記多項式パラメータから構成される当てはめ曲線と、原時系列データの対応部分との差演算により求められる残差時系列データを次の下位レイヤーのセグメントデータとして設定する処理を行う2系列の時系列データの解析装置を提案する。 In the present invention, in the above analysis device, the analysis of the segment data in the data analysis means is performed sequentially from the highest layer to the lowest layer. At this time, the segment data of the upper layer segment is In addition, in order to prevent the number of data points per segment from exceeding the preset maximum number of data points, the processing is performed to reduce the number of data points of the time series data to be analyzed by bunching processing and configure segment data, A temporary storage area is configured corresponding to each segment of the layer, and the polynomial parameter obtained for the upper segment is stored in the temporary storage area of each segment of each lower layer, and is configured from the above polynomial parameter. By calculating the difference between the fitting curve and the corresponding part of the original time series data. The residual time series data that is to propose an analysis device of the time-series data of the two series which performs a process of setting as the segment data of the next lower layer.
また本発明では、上記の解析装置において、データ解析手段により求めた解析結果を評価し、条件に適合する場合に記憶手段に記憶する解析結果評価手段と、上記パラメータから得られるパワースペクトル密度と最大エントロピー法を用いて求める残差のパワースペクトル密度を合成して時系列データの再構成スペクトルを求める再構成スペクトル計算手段と、この計算手段を利用して、上記記憶手段から取り出した上記パラメータと残差により、解析対象の時系列データに対応する一般化三角多項式を再構成する解析結果構成手段とを備えた2系列の時系列データの解析装置を提案する。 In the present invention, in the above analysis apparatus, the analysis result obtained by the data analysis means is evaluated, and when the condition is met, the analysis result evaluation means for storing in the storage means, the power spectral density obtained from the parameters and the maximum A reconstructed spectrum calculating means for obtaining a reconstructed spectrum of time series data by synthesizing the power spectrum density of the residual obtained using the entropy method, and using the calculating means, the parameters extracted from the storing means and the remaining A two-series time-series data analysis device is provided, which includes an analysis result construction means for reconstructing a generalized trigonometric polynomial corresponding to the time-series data to be analyzed based on the difference.
また本発明では、上記の解析装置において、解析結果構成手段は、全ての記憶領域に記憶されている基本セグメントのパラメータと残差とから全データ領域を一括して再構成スペクトルを求めて記憶する一括解析結果構成手段と、各記憶領域に記憶されている基本セグメントのパラメータと残差とから、基本セグメント毎に再構成スペクトルを求めて記憶するセグメント解析結果構成手段とから構成すると共に、基本セグメント毎に頂位位相の周波数依存性を解析するための頂位位相推移計算手段を構成し、この頂位位相推移計算手段により計算した各基本セグメント毎の頂位位相の周波数依存をセグメント解析結果構成手段による再構成スペクトルと共に、基本セグメントの代表時刻の属性として記憶する2系列の時系列データの解析装置を提案する。 Also, in the present invention, in the above analysis apparatus, the analysis result configuration means obtains and stores a reconstructed spectrum for all data areas collectively from the parameters and residuals of the basic segments stored in all the storage areas. The basic analysis segment is composed of a batch analysis result configuration unit and a segment analysis result configuration unit that obtains and stores a reconstruction spectrum for each basic segment from the parameters and residuals of the basic segment stored in each storage area. Configure the top phase transition calculation means to analyze the frequency dependence of the top phase every time, and configure the segment analysis result for the frequency dependence of the top phase for each basic segment calculated by this top phase transition calculation means Apparatus for analyzing two series of time series data stored as an attribute of a representative time of a basic segment together with a reconstructed spectrum by means Proposed.
また本発明では、上記の解析装置において、再構成スペクトル計算手段は、得られたパラメータからパワースペクトル密度を求める計算において、一般化三角多項式の夫々の項の夫々につき、δ函数の原型を用いて全周波数帯に裾の広がる夫々一つのスペクトルピークを構成することにより、パラメータから得られる各周波数に対応する孤立ピーク群を構成し、その重畳として時系列データの再構成スペクトルを求める構成とする2系列の時系列データの解析装置を提案する。 Further, in the present invention, in the above analysis device, the reconstructed spectrum calculating means uses the prototype of the δ function for each of the terms of the generalized triangular polynomial in the calculation for obtaining the power spectral density from the obtained parameters. By configuring each spectrum peak having a wide base in the entire frequency band, an isolated peak group corresponding to each frequency obtained from the parameters is formed, and a reconstructed spectrum of time series data is obtained as a superposition thereof. A time series data analysis device is proposed.
また本発明では、上記の構成において、δ函数の原型として、余弦函数に微弱な雑音を加えた時系列をラグ2にて最大エントロピー法のスペクトル密度の計算を行って得られたスペクトルを利用する2系列の時系列データの解析装置を提案する。
In the present invention, in the above configuration, as a prototype of the δ function, a spectrum obtained by calculating the spectral density of the maximum entropy method at
また本発明では、解析対象の2系列の時系列データを、夫々解析して、その一般化三角多項式パラメータを求める単一系列解析手順と、単一系列解析手順により求められた2系列の時系列データの一般化三角多項式パラメータから、それらの時系列データの間の相互相関を定量化する諸量を求める相互相関解析手順とから構成し、単一系列解析手順は、解析対象の時系列データを解析するシーケンスを決定する解析シーケンス決定手順と、解析シーケンスに従い、最大エントロピー法及び非線形最小自乗法によりデータを解析して一般化三角多項式のパラメータと残差を求めるデータ解析手順とから構成し、相互相関解析手順は、一般化三角多項式パラメータからパワースペクトル密度を求めるパワースペクトル密度計算手順と、クロススペクトル(コスペクトル、クオドスペクトル)とコヒーレンスを含む相互相関を定量化する諸量を求める相互相関諸量計算手順と、相互相関諸量計算手順により求められたコヒーレンスを理論値と比較して、それらの差が所定値以下の場合に、その差をなくすようにクロススペクトルを補正する整合性評価補正手順と、整合性評価補正手順において評価した差が所定値以上の場合に、単一系列解析手順に再解析の指令を出力し、所定値以下の場合には、相互相関諸量計算手順により計算された諸量と、パワースペクトル密度計算手順により計算されたパワースペクトル密度及び整合性評価補正手順により補正されたクロススペクトルを出力する分岐出力手順を構成した2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。 In the present invention, two series of time series data to be analyzed are analyzed, and a single series analysis procedure for obtaining the generalized trigonometric polynomial parameters, and two series of time series obtained by the single series analysis procedure. The cross-correlation analysis procedure is used to obtain various quantities that quantify the cross-correlation between the time series data from the generalized trigonometric polynomial parameters of the data. The analysis sequence determination procedure for determining the sequence to be analyzed, and the data analysis procedure for analyzing the data by the maximum entropy method and the nonlinear least square method according to the analysis sequence to obtain the parameters and residuals of the generalized triangular polynomial, The correlation analysis procedure includes a power spectrum density calculation procedure for obtaining the power spectrum density from the generalized trigonometric polynomial parameters, and a cross spectrum. Cross-correlation quantities calculation procedure for obtaining quantities that quantify cross-correlation including Torr (cospectrum, quad spectrum) and coherence, and comparing the coherence obtained by the cross-correlation quantity calculation procedure with theoretical values. If the difference is less than or equal to a predetermined value, the consistency evaluation correction procedure for correcting the cross spectrum so as to eliminate the difference, and the single series analysis procedure when the difference evaluated in the consistency evaluation correction procedure is greater than or equal to the predetermined value If the re-analysis command is output to the value below the predetermined value, the amount calculated by the cross-correlation amount calculation procedure and the power spectrum density and consistency evaluation correction procedure calculated by the power spectrum density calculation procedure are used. Computer-readable recording of two series of time series data analysis programs that constitute a branch output procedure that outputs a corrected cross spectrum Suggest Do recording medium.
また本発明では、上記の記録媒体において、パワースペクトル密度計算手順は、一般化三角多項式パラメータからパワースペクトル密度を求める計算において、一般化三角多項式の夫々の項の夫々につき、δ函数の原型を用いて全周波数帯に裾の広がる夫々一つのスペクトルピークを構成することにより、パラメータから得られる各周波数に対応する孤立ピーク群を構成し、その重畳として時系列データのスペクトルを求める構成とする2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。 In the present invention, in the above recording medium, the power spectral density calculation procedure uses the prototype of the δ function for each of the terms of the generalized triangular polynomial in the calculation for obtaining the power spectral density from the generalized triangular polynomial parameter. By configuring each spectrum peak that spreads in the entire frequency band, an isolated peak group corresponding to each frequency obtained from the parameters is configured, and two sequences are configured to obtain the spectrum of the time series data as its superposition We propose a computer-readable recording medium in which a time series data analysis program is recorded.
また本発明では、上記の記録媒体において、δ函数の原型として、余弦函数に微弱な雑音を加えた時系列をラグ2にて最大エントロピー法のスペクトル密度の計算を行って得られたスペクトルを利用する2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。
In the present invention, in the above recording medium, a spectrum obtained by calculating the spectral density of the maximum entropy method with
また本発明では、上記の記録媒体において、解析シーケンス決定手順は、解析対象の時系列データを、予め設定された最大データ点数を越えないようにセグメントに分割する分割処理過程と、各セグメントに対してデータ構造の均一性を評価する評価処理過程と、この評価に基づき、不均一のセグメントデータを均一の複数のセグメントに分割する分割処理過程と、分割された均一のセグメントの夫々を基本セグメントとして、それらに対して、記憶手段に記憶領域を設定する記憶領域設定処理過程を備えている2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。 In the present invention, in the above recording medium, the analysis sequence determination procedure includes a division process for dividing the time-series data to be analyzed into segments so as not to exceed a preset maximum number of data points, and for each segment. The evaluation process for evaluating the uniformity of the data structure, the division process for dividing the non-uniform segment data into a plurality of uniform segments based on this evaluation, and each of the divided uniform segments as a basic segment Accordingly, a computer-readable recording medium recording a two-series time-series data analysis program having a storage area setting process for setting a storage area in the storage means is proposed.
また本発明では、上記の記録媒体において、セグメントは、解析対象の時系列データの全期間データに対応する最上位レイヤーのセグメントから、順次2分割されて構成される、下位のレイヤーの夫々に含まれる2の累乗のセグメントと、最下位のレイヤーにおいて、データ構造が変化するセグメントを、データ構造の変化している個所において分割して順次形成されるサブレイヤーの各セグメントとから構成される2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。 According to the present invention, in the recording medium, the segment is included in each of the lower layers configured by sequentially dividing the segment into the uppermost layer segment corresponding to the entire period data of the time-series data to be analyzed. 2 series comprising segments of powers of 2 and sub-layer segments sequentially formed by dividing a segment whose data structure changes in the lowest layer at a location where the data structure changes We propose a computer-readable recording medium in which a time series data analysis program is recorded.
また本発明では、上記の記録媒体において、データ構造の均一性は、下式による最大ラグ値までの最大エントロピー法のスペクトル密度の計算におけるバーグの係数Pmの、ラグ値に対する推移により評価する2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。
また本発明では、上記の記録媒体において、解析シーケンス決定手順には、不等間隔の時系列データを等間隔化処理する等間隔化処理過程を設けた2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。 Further, in the present invention, in the above recording medium, the analysis sequence determination procedure records a two-sequence time-series data analysis program provided with an equal-interval processing process for equal-interval processing of unequal-interval time-series data. A computer-readable recording medium is proposed.
また本発明では、上記の記録媒体において、データ解析手順は、解析対象の時系列データから取り出したセグメントデータに対して下式による最大ラグ値までの最大エントロピー法のスペクトル密度の計算におけるバーグの係数Pmの、ラグ値に対する推移を求める計算過程と、求めた推移から最大エントロピー法によるパワースペクトル密度の計算に用いるラグ値を決定するラグ値決定過程と、決定されたラグ値により上記パワースペクトル密度を求めるパワースペクトル密度計算過程と、求められたパワースペクトル密度から主要ピークの周波数と、そのパワーを抽出する抽出過程と、抽出された周波数とパワーとから非線形最小自乗法により当てはめ曲線のパラメータを求めるパラメータ導出過程と、求められた当てはめ曲線とセグメントデータとを比較して一致度を評価し、所定の一致度以上の場合には上記当てはめ曲線のパラメータをセグメントに対応する記憶手段に記憶して次の過程に移行し、所定の一致度に達しない場合には上記当てはめ曲線のパラメータを破棄して上記ラグ値を決定する過程に移行する評価過程とを備えており、ラグ値決定過程は、設定される解析の精度に応じて、より高精度の解析において、より大きいラグ値を選択するように構成し、抽出過程は、設定される解析の精度に応じて、より高精度の解析において、より多くのピークを主要ピークとして抽出するように構成すると共に、評価過程は、設定される解析の精度に応じて、より高精度の解析において、より高い一致度となるように構成し、相互相関解析手順からの再解析指令により、高精度の解析が設定されるように構成した2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。
また本発明では、上記の記録媒体において、データ解析手順におけるセグメントデータの解析は、最上位のレイヤーから最下位のレイヤーに向かって順次行うものとし、この際、上位のレイヤーのセグメントのセグメントデータは、セグメント当たりのデータ点数が予め設定された最大データ点数を越えないように、解析対象の時系列データのデータ点数をバンチング処理により低減してセグメントデータを構成する処理を行う構成とすると共に、各レイヤーの各セグメントに対応して一時記憶領域を構成し、上位側のセグメントに対して求められた多項式パラメータを下位の各レイヤーの各セグメントの一時記憶領域に記憶すると共に、上記多項式パラメータから構成される当てはめ曲線と、原時系列データの対応部分との差演算により求められる残差時系列データを次の下位レイヤーのセグメントデータとして設定する処理を行う2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。 In the present invention, in the recording medium, the segment data analysis in the data analysis procedure is performed sequentially from the highest layer to the lowest layer. At this time, the segment data of the upper layer segment is In addition, in order to prevent the number of data points per segment from exceeding the preset maximum number of data points, the processing is performed to reduce the number of data points of the time series data to be analyzed by bunching processing and configure segment data, A temporary storage area is configured corresponding to each segment of the layer, and the polynomial parameter obtained for the upper segment is stored in the temporary storage area of each segment of each lower layer, and is configured from the above polynomial parameter. By calculating the difference between the fitting curve and the corresponding part of the original time series data. The residual time series data that is to propose a computer-readable recording medium storing an analysis program of the time-series data of the two series which performs a process of setting as the segment data of the next lower layer.
また本発明では、上記の記録媒体において、データ解析手順により求めた解析結果を評価し、条件に適合する場合に記憶手段に記憶する解析結果評価手順と、上記パラメータから得られるパワースペクトル密度と最大エントロピー法を用いて求める残差のパワースペクトル密度を合成して時系列データの再構成スペクトルを求める再構成スペクトル計算手順と、この計算手順を利用して、上記記憶手段から取り出した上記パラメータと残差により、解析対象の時系列データに対応する一般化三角多項式を再構成する解析結果構成手順とを備えた2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。 In the present invention, in the recording medium, the analysis result obtained by the data analysis procedure is evaluated, and the analysis result evaluation procedure stored in the storage means when the condition is satisfied, the power spectral density obtained from the parameters, and the maximum A reconstructed spectrum calculation procedure for obtaining a reconstructed spectrum of time series data by synthesizing the power spectrum density of the residual obtained using the entropy method, and using the calculation procedure, the parameters extracted from the storage means and the residual According to the difference, a computer-readable recording medium recording a two-series time series data analysis program including an analysis result configuration procedure for reconstructing a generalized trigonometric polynomial corresponding to time series data to be analyzed is proposed.
また本発明では、上記の記録媒体において、解析結果構成手順は、全ての記憶領域に記憶されている基本セグメントのパラメータと残差とから全データ領域を一括して再構成スペクトルを求めて記憶する一括解析結果構成手順と、各記憶領域に記憶されている基本セグメントのパラメータと残差とから、基本セグメント毎に再構成スペクトルを求めて記憶するセグメント解析結果構成手順とから構成すると共に、基本セグメント毎に頂位位相の周波数依存性を解析するための頂位位相推移計算手順を構成し、この頂位位相推移計算手順により計算した各基本セグメント毎の頂位位相の周波数依存をセグメント解析結果構成手順による再構成スペクトルと共に、基本セグメントの代表時刻の属性として記憶する2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。 In the present invention, in the above recording medium, the analysis result composition procedure obtains and stores the reconstructed spectrum for all data areas collectively from the basic segment parameters and residuals stored in all the storage areas. It consists of a batch analysis result configuration procedure and a segment analysis result configuration procedure for obtaining and storing a reconstructed spectrum for each basic segment from the parameters and residuals of the basic segment stored in each storage area. Configure the top phase transition calculation procedure to analyze the frequency dependence of the top phase every time, and configure the segment analysis result for the frequency dependence of the top phase for each basic segment calculated by this top phase transition calculation procedure An analysis program for two series of time series data that is stored as an attribute of the representative time of the basic segment along with the reconstructed spectrum by the procedure It proposes a computer-readable recording medium recording a ram.
また本発明では、上記の記録媒体において、再構成スペクトル計算手順は、得られたパラメータからパワースペクトル密度を求める計算において、一般化三角多項式の夫々の項の夫々につき、δ函数の原型を用いて全周波数帯に裾の広がる夫々一つのスペクトルピークを構成することにより、パラメータから得られる各周波数に対応する孤立ピーク群を構成し、その重畳として時系列データの再構成スペクトルを求める構成とする2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。 In the present invention, in the recording medium described above, the reconstruction spectrum calculation procedure uses the prototype of the δ function for each of the terms of the generalized trigonometric polynomial in the calculation for obtaining the power spectral density from the obtained parameters. By configuring each spectrum peak having a wide base in the entire frequency band, an isolated peak group corresponding to each frequency obtained from the parameters is formed, and a reconstructed spectrum of time series data is obtained as a superposition thereof. A computer-readable recording medium on which a time series data analysis program is recorded is proposed.
また本発明では、上記の記録媒体において、δ函数の原型として、余弦函数に微弱な雑音を加えた時系列をラグ2にて最大エントロピー法のスペクトル密度の計算を行って得られたスペクトルを利用する2系列の時系列データの解析プログラムを記録したコンピュータ読み取り可能な記録媒体を提案する。
In the present invention, in the above recording medium, a spectrum obtained by calculating the spectral density of the maximum entropy method with
請求項1、14の発明によれば、2系列の時系列データの夫々を単一系列解析手段において、まず最大エントロピー法を適用して高精度なパワースペクトル密度を求めると共に、その顕著なピークから、一般化三角多項式を解く初期値としての、項数と、各三角項の周期値を求め、次いで、これらの初期値により非線形最小自乗法を用いて、各三角項のパラメータを求めると共に、残差時系列データを求め、これらを加算して時系列データを再構成するという解析を行うことにより、何らの数値モデルも想定せず、夫々の時系列データのみを用いて高精度の解析を行え、これらの高精度の解析によって得られた2系列の時系列データの一般化三角多項式パラメータにより、相互相関解析手段において、2系列の時系列データの間の相互相関を定量化する諸量を求めることができる。
According to the inventions of
特に上述した発明によれば、相互相関解析手段においては、一般化三角多項式パラメータにより求められたコヒーレンスを理論値と比較して整合性を評価し、所定以上の整合性が得られない場合には、単一系列解析手段に対して再解析指令を発して、より高精度の解析を行わせるという、フィードバックと繰り返し動作を有するため、2系列の時系列データの間の相互相関を定量化する諸量を非常に高精度に求めることができる。 In particular, according to the above-described invention, the cross-correlation analysis means evaluates the consistency by comparing the coherence obtained by the generalized trigonometric polynomial parameter with the theoretical value, and when the consistency exceeding a predetermined value cannot be obtained. Because it has a feedback and repetitive action to issue a reanalysis command to a single series analysis means to perform a more accurate analysis, various measures to quantify the cross-correlation between two series of time series data The quantity can be determined with very high accuracy.
請求項2、3、12、13、15、16、25及び26の発明によれば、理論形式では有限離散値を扱うコンピュータで構成することができない線スペクトル列となる時系列データのスペクトルを、ピーク周波数から離れるに従って急減する連続関数、例えば余弦函数に微弱な雑音を加えた時系列をラグ2にて最大エントロピー法のスペクトル密度の計算を行って得られたスペクトルをδ函数の原型として利用して再構成することにより、コンピュータを利用して、理論と整合する時系列データの解析を行うことができる。また、これにより、2つの時系列データのクロススペクトルの絶対値と等しくなる、2つの時系列データの積を構成することも可能となる。
According to the inventions of
請求項4、17の発明によれば、解析シーケンス決定手段又は手順により、解析対象の時系列データを、その全体のデータ点数に応じて分割して解析装置の演算能力に適合させると共に、データ構造が変化している部分においてセグメントに分割することにより、セグメント中にデータ構造が変化している部分を無くすことにより、これらの均一なセグメントにつき、データ解析手段又は手順により高精度に解析を行うことができる。
According to the inventions of
請求項5、18の発明によれば、セグメントは、解析シーケンス決定手段又は手順により、解析対象の時系列データの全期間データに対応する、最上位レイヤーのセグメントから、順次2分割されて構成される、下位のレイヤーの夫々に含まれる2の累乗のセグメントと、最下位のレイヤーにおいて、データ構造が変化するセグメントを、データ構造の変化している個所において分割して順次形成されるサブレイヤーの各セグメントとから構成するため、解析装置の演算能力と時系列データのデータ構造の両方に対応して、合理的に時系列データのセグメント化を行うことができ、自動化も容易である。
According to the inventions of
請求項6、19の発明によれば、データ構造の均一性を容易に評価することができる。
According to the inventions of
請求項7、20の発明によれば、原時系列データが不等間隔の場合にも、等間隔化することにより、最大エントロピー法のスペクトル密度の計算を行うことができる。一方、非線形最小自乗法により当てはめ曲線のパラメータを求めるパラメータの導出においては、必要に応じて、等間隔化した時系列データか、不等間隔の原時系列データをそのまま使用することができる。 According to the seventh and twentieth aspects of the present invention, even when the original time series data is unequal, the spectral density of the maximum entropy method can be calculated by equalizing the original time series data. On the other hand, in the derivation of the parameters for obtaining the parameters of the fitting curve by the non-linear least square method, it is possible to use time-series data with equal intervals or original time-series data with unequal intervals as necessary.
請求項8、21の発明によれば、データ解析手段又は手順において、求められた当てはめ曲線とセグメントデータとを比較して一致度を評価し、所定の一致度以上の場合には上記当てはめ曲線のパラメータをセグメントに対応する記憶手段に記憶して次の過程に移行し、所定の一致度に達しない場合には上記当てはめ曲線のパラメータを破棄して上記ラグ値を決定する過程に移行する評価過程とを備えており、パラメータの導出と、その評価を繰り返すことにより、単一系列解析手段自体において、より高精度の当てはめ曲線を求めることができる。
According to the inventions of
更にこれらの発明によれば、ラグ値決定過程は、設定される解析の精度に応じて、より高精度の解析において、より大きいラグ値を選択するように構成され、また抽出過程は、設定される解析の精度に応じて、より高精度の解析において、より多くのピークを主要ピークとして抽出するように構成され、更に評価過程は、設定される解析の精度に応じて、より高精度の解析において、より高い一致度となるように構成されているので、相互相関解析手段からの再解析指令毎に、より高精度の解析が行われる。勿論、解析の精度は、相互相関解析手段からの再解析指令のみでなく、操作者が設定可能とすることもできる。 Furthermore, according to these inventions, the lag value determination process is configured to select a larger lag value in a higher-accuracy analysis according to the accuracy of the set analysis, and the extraction process is set. Depending on the accuracy of analysis, it is configured to extract more peaks as main peaks in higher accuracy analysis, and the evaluation process is more accurate depending on the set accuracy of analysis. In FIG. 1, the higher the degree of coincidence, the more accurate analysis is performed for each re-analysis command from the cross-correlation analysis means. Of course, the accuracy of the analysis can be set not only by the re-analysis command from the cross-correlation analysis means but also by the operator.
請求項9、22の発明によれば、解析対象の時系列データの全体のデータ点数を一括して処理するのには、解析装置の演算能力が足りない場合であっても、バンチング処理によりデータ点数を低減することにより、最上位のレイヤーのセグメントデータ、即ち、解析対象の時系列データの全体の一括した解析が可能であり、従って時系列データの全体に渡る長周期の解析が可能であると共に、上位のレイヤーの解析により得られた当てはめ曲線と、原時系列データの対応部分との差演算により求められる残差時系列データを次の下位レイヤーのセグメントデータとして設定して、順次解析を行うので、最下位又はサブレイヤーのセグメントデータにより、バンチング処理によるデータ点数の低減がない、そのままの時系列データに対しての解析を行うことができる。 According to the ninth and twenty-second aspects of the present invention, the entire data points of the time-series data to be analyzed can be collectively processed even if the analysis apparatus has insufficient computing power, and the data is obtained by bunching processing. By reducing the number of points, the segment data of the highest layer, that is, the entire time-series data to be analyzed can be collectively analyzed, and therefore the long-period analysis over the entire time-series data is possible. At the same time, set the residual time series data obtained by calculating the difference between the fitting curve obtained by the analysis of the upper layer and the corresponding part of the original time series data as the segment data of the next lower layer, and analyze sequentially. As a result, the solution for the time-series data as it is without the reduction in the number of data points due to the bunching process due to the segment data of the lowest or sublayer. It can be carried out.
請求項10、23の発明によれば、解析結果評価手段又は手順において、得られた当てはめ曲線が、解析対象の時系列データの対応範囲と、所定の一致度に達しない場合にはデータを破棄して、セグメントデータを次のレイヤーのセグメントに移行して、処理に供されるので、上位のレイヤーにおいて一致度の高い当てはめ曲線が求まらない場合でも、不都合が生じず、下位のレイヤーにおいての処理により、一致度の高い当てはめ曲線を得ることができる。
According to the inventions of
請求項11、24の発明によれば、単一系列解析手段において、解析結果構成手段又は手順により、解析対象の時系列データを、一括した解析結果と、各セグメント毎の解析結果として出力することができる。
According to the inventions of
次に本発明に係る2系列の時系列データの解析装置の最良の形態を、添付図面を参照して説明する。
まず、図1に示すように本発明に係る統合構造解析手段1は、解析対象の2系列の時系列データを、夫々解析して、その一般化三角多項式パラメータを求める2系統の単一系列解析手段2a,2bと、夫々の単一系列解析手段2a,2bにより求められた2系列の時系列データの一般化三角多項式パラメータから、2系列の時系列データの間の相互相関を定量化する諸量を求める相互相関解析手段3とから構成されている。夫々の単一系列解析手段2a,2bにより求められた2系列の時系列データの一般化三角多項式パラメータは、図中実線矢印に示すように相互相関解析手段3に入力され、一方、図中破線矢印で示すように相互相関解析手段3から夫々の単一系列解析手段2a,2bに必要に応じて再解析指令が発せられる。尚、図1において、符号4a,4bは2系列の夫々の時系列データ、符号5a,5b,5cは本発明の解析装置から出力される諸量を示すものである。
Next, the best mode of a two-series time-series data analyzing apparatus according to the present invention will be described with reference to the accompanying drawings.
First, as shown in FIG. 1, the integrated structural analysis means 1 according to the present invention analyzes two series of time series data to be analyzed, and obtains a generalized trigonometric polynomial parameter for two systems of single series analysis. Various means for quantifying the cross-correlation between the two series of time series data from the generalized trigonometric polynomial parameters of the two series of time series data obtained by the
単一系列解析手段2a,2bは、上述した解析手法を適用するため、解析対象の時系列データ(以降、必要に応じて原時系列データと称する)を一般化三角多項式により表現して解析を行うように構成されており、概ね、原時系列データに対して最大エントロピー法を適用して、高精度なパワースペクトル密度(PSD)を求め、その顕著なピークから、一般化三角多項式を解く初期値としての、項数と、各三角項の周期値を求める過程と、これらの初期値により、非線形最小自乗法を用いて、各三角項のパラメータを求めると共に、残差時系列データを求め、これらを加算して時系列データを再構成する過程を実行するように構成されている。 In order to apply the above-described analysis method, the single series analysis means 2a and 2b express the time series data to be analyzed (hereinafter referred to as original time series data as necessary) by a generalized trigonometric polynomial. In general, the maximum entropy method is applied to the original time series data to obtain high-accuracy power spectral density (PSD), and the initial stage of solving the generalized trigonometric polynomial from its prominent peaks The process of obtaining the number of terms as a value and the period value of each triangular term and the initial values thereof, using the nonlinear least squares method, obtain the parameters of each triangular term and obtain the residual time series data, A process of adding these to reconstruct time-series data is executed.
即ち、本発明に係る単一系列解析手段2a,2bにおいては、原時系列データは、次式の一般化三角多項式で表現される。
図2は単一系列解析手段2a,2bの実施の形態を模式的に表した全体構成の説明図であり、図中矢印により処理の流れを模式的に表している。
図2に示すように、この実施の形態に係る単一系列解析手段2a,2bは、原時系列データを解析するシーケンスを決定する解析シーケンス決定手段6と、この解析シーケンス決定手段6により決定された解析シーケンスに従い、最大エントロピー法及び非線形最小自乗法によりデータを解析して一般化三角多項式のパラメータと残差を求めるデータ解析手段7と、データ解析手段7により求めた解析結果を評価し、条件に適合する場合に、予め設定された記憶手段に記憶する解析結果評価手段8と、これらの解析シーケンス決定手段6と、データ解析手段7と、解析結果評価手段8を統合的に制御して、必要な繰り返し演算等を行わせるための制御手段9と、最大エントロピー法を用いて一般化三角多項式の再構成スペクトルを求める再構成スペクトル計算手段10と、頂位位相の推移を計算する頂位位相推移計算手段11と、出力するデータの準備等に係る計算を行う、その他計算手段12と、これらの各計算手段10〜12を用いて、出力するデータ解析手段7の解析結果を構成する解析結果構成手段13a,13bとから構成されている。
FIG. 2 is an explanatory diagram of the overall configuration schematically showing the embodiment of the single series analyzing means 2a, 2b, and the flow of processing is schematically shown by arrows in the figure.
As shown in FIG. 2, the single series analysis means 2a, 2b according to this embodiment is determined by an analysis sequence determination means 6 for determining a sequence for analyzing original time series data and the analysis sequence determination means 6. In accordance with the analysis sequence described above, the data is analyzed by the maximum entropy method and the nonlinear least square method to obtain the parameters and residuals of the generalized trigonometric polynomial, the analysis results obtained by the data analysis means 7 are evaluated, , The analysis result evaluation means 8 stored in a preset storage means, the analysis sequence determination means 6, the data analysis means 7, and the analysis result evaluation means 8 are integratedly controlled, A control means 9 for performing necessary iterative calculations and the like, and a reconstruction for obtaining a reconstructed spectrum of a generalized triangular polynomial using the maximum entropy method. The spectrum calculation means 10, the top phase transition calculation means 11 for calculating the transition of the top phase, the other calculation means 12 for performing calculations related to the preparation of data to be output, and the calculation means 10 to 12 And the analysis result constituting means 13a and 13b constituting the analysis result of the data analyzing means 7 to be output.
後述するように、解析シーケンス決定手段6は、解析シーケンスの一つとして、原時系列データの全期間データのデータ点数やデータ構造等に応じて、それを複数レイヤーに構成した、データ構造の均一な複数のセグメントに分割すると共に、分割された均一なセグメントの夫々を基本セグメントとして、それらに対応して記憶手段に記憶領域を設定する過程を有し、この分割されたセグメントに取り込まれた原時系列データについてデータ解析手段7において解析されるものである。そしてデータ解析手段7において各セグメントの原時系列データについて得られた解析結果は、対象とする時系列データ全体を一括する一括解析結果と、セグメント解析結果として構成されるものであり、このために上記解析結果構成手段13a,13bは、前者に対応する一括解析結果構成手段13aと、後者に対応するセグメント解析結果構成手段13bとから構成されている。 As will be described later, the analysis sequence determination means 6 has a uniform data structure in which the analysis sequence determination means 6 is configured as a plurality of layers according to the data points, data structures, etc. of the whole period data of the original time series data as one of the analysis sequences. And dividing each of the divided uniform segments into basic segments, and setting a storage area in the storage means corresponding to them, and the original captured in the divided segments. The time series data is analyzed by the data analysis means 7. The analysis results obtained for the original time-series data of each segment in the data analysis means 7 are configured as a batch analysis result for collecting the entire target time-series data and a segment analysis result. The analysis result configuration means 13a and 13b are composed of collective analysis result configuration means 13a corresponding to the former and segment analysis result configuration means 13b corresponding to the latter.
図3はこの実施の形態の単一系列解析手段2a,2bから出力される解析結果を模式的に示すものであり、原時系列データの全期間に渡ってデータ構造が保たれる一括の解析結果14と、原時系列データの期間内でデータ構造が時々刻々と変動し、原時系列データの全期間を網羅する複数のセグメント毎に得られる解析結果15とから構成される。従って、前者の解析結果14は時間依存がなく、後者の解析結果15は時間依存を有するものである。
FIG. 3 schematically shows the analysis results output from the single series analysis means 2a, 2b of this embodiment, and a batch analysis in which the data structure is maintained over the entire period of the original time series data. The
解析結果14において、一般化三角多項式パラメータaは、上述した一般化三角多項式の項数、振幅、周期及び頂位位相であり、これらのパラメータaから、周期的に変動する曲線、即ち、多項式パラメータによる最適あてはめ曲線bと、そのスペクトル、即ち、多項式パラメータによる再構成スペクトルcが構成され、これに残差と残差のスペクトルdによる残差成分を加算することにより、原時系列データとそのスペクトルが再構成される。
In the
一方、解析結果15においては、一般化三角多項式パラメータeと、多項式パラメータによる再構成スペクトルfと、多項式パラメータによる頂位位相gとが夫々各セグメント毎に構成され、従って再構成スペクトルfと頂位位相gの時間変動が構成される。
On the other hand, in the
次に、単一系列解析手段2a,2bの、上述した各構成要素の動作を説明する。まず図4、図5は解析シーケンス決定手段6における処理の流れを示す流れ図であり、図6は解析シーケンス決定手段6により決定されたレイヤーとセグメントの構成例を模式的に示すものである。 Next, the operation of each component described above of the single series analysis means 2a and 2b will be described. 4 and 5 are flowcharts showing the flow of processing in the analysis sequence determination means 6, and FIG. 6 schematically shows a configuration example of layers and segments determined by the analysis sequence determination means 6.
まず、単一系列解析手段2a,2bにステップS1において原時系列データが入力されると、解析シーケンス決定手段6は、まずステップS2において、一つのセグメントが配置されるレイヤーを設定する。この一つのセグメントが配置されるレイヤーは必ず設定されるものであるので、その設定時点はステップS2の時点でなくとも良い。 First, when original time series data is input to the single series analysis means 2a and 2b in step S1, the analysis sequence determination means 6 first sets a layer in which one segment is arranged in step S2. Since the layer in which this one segment is arranged is always set, the setting time may not be the time of step S2.
次にステップS3では、入力された原時系列データの全期間のデータ点数が、予め設定されている最大データ点数を越えているか、否かを判定する。
この最大データ点数は、後述する基本セグメントに取り込んで解析が行われる最大データ点数であるが、一般に、各セグメントに取り込んで同時に解析するデータ点数が小さすぎるとデータの構造の抽出が困難となり、逆に大きすぎると過剰なデータ構造が抽出されてしまうため、最適なデータ点数とするのが高精度の解析結果を得るためには好ましい。一方、解析に最適なデータ点数は、取り込んだデータの周期性が強いか否かによっても異なり、また単一系列解析手段2a,2bを構成する演算手段の演算速度や演算精度等によってもデータ点数が制限される。
Next, in step S3, it is determined whether or not the total number of data points of the input original time series data exceeds the preset maximum number of data points.
This maximum number of data points is the maximum number of data points that can be analyzed by importing into the basic segment, which will be described later. Generally, if the number of data points that are imported into each segment and analyzed at the same time is too small, it will be difficult to extract the structure of the data. If it is too large, an excessive data structure is extracted. Therefore, it is preferable to obtain an optimal number of data points in order to obtain a highly accurate analysis result. On the other hand, the optimal number of data points for analysis varies depending on whether or not the captured data has a strong periodicity, and the number of data points also depends on the calculation speed and calculation accuracy of the calculation means constituting the single series analysis means 2a and 2b. Is limited.
このようなことから、最大データ点数としては、それまでの多数の時系列データの解析により蓄積された知見等から得られた値が設定される。尚、この最大データ点数は、後述する各評価手段による評価の結果をフィードバックして変更可能とし、変更された最大データ点数に基づいてステップS3及びそれ以降の処理を行うようにすることもできる。また最大データ点数は、操作者により設定可能とすることもできる。 For this reason, as the maximum number of data points, a value obtained from knowledge accumulated by analysis of a large number of time series data so far is set. Note that this maximum data score can be changed by feeding back the results of evaluation by each evaluation means to be described later, and step S3 and subsequent processing can be performed based on the changed maximum data score. The maximum number of data points can be set by the operator.
ステップS3において、入力された時系列データの全体のデータ点数が、予め設定されている最大データ点数を越えていないと判定された場合には、一つのセグメントが配置された一つのレイヤーが設定された状態で、結合子Aを経て図5のステップS6に移行する。 If it is determined in step S3 that the total number of data points of the input time-series data does not exceed the preset maximum number of data points, one layer in which one segment is arranged is set. In this state, the process proceeds to step S6 in FIG.
一方、ステップS3において、入力された時系列データの全体のデータ点数が、予め設定されている最大データ点数を越えていると判定された場合には、ステップS4に移行し、そのステップS4において、それまでに設定された上位のレイヤーの下位にレイヤーを追加し、その下位のレイヤーには、上位のレイヤーのセグメント数の2倍のセグメントを配置して、ステップS5に移行する。 On the other hand, if it is determined in step S3 that the total number of data points of the input time-series data exceeds the preset maximum number of data points, the process proceeds to step S4. In step S4, A layer is added to the lower layer of the upper layer set so far, and a segment twice as many as the number of segments of the upper layer is arranged in the lower layer, and the process proceeds to step S5.
ステップS5では、入力された原時系列データの全期間データが、ステップS4において配置された全セグメントに等分されて取り込まれた場合の各セグメント毎のデータ点数が上記最大データ点数を越えているか、否かを判定する。 In step S5, whether the number of data points for each segment when the whole period data of the input original time series data is equally divided into all the segments arranged in step S4 exceeds the maximum number of data points. Determine whether or not.
ステップS5において、各セグメントのデータ点数が、上記最大データ点数を越えていると判定された場合には、ステップS4に移行して、更に下位のレイヤーが追加され、一方、各セグメントのデータ点数が、上記最大データ点数を越えていないと判定された場合には、結合子Aを経て図5のステップS6に移行する。 If it is determined in step S5 that the number of data points of each segment exceeds the maximum data point, the process proceeds to step S4, where a lower layer is added, while the number of data points of each segment is If it is determined that the maximum number of data points has not been exceeded, the process proceeds to step S6 in FIG.
以上の各ステップにおける処理を経て、原時系列データを、予め設定された最大データ点数を越えないようにセグメントに取り込むための、レイヤーとセグメントの構成が得られ、構成された単数又は順次配列された複数のセグメントが、原時系列データの全期間を網羅する。従って以上のステップの処理は、上述した分割処理過程に相当する。 Through the processing in each of the above steps, the structure of layers and segments for obtaining original time series data into a segment so as not to exceed the preset maximum number of data points is obtained, and the configured single or sequential Multiple segments cover the entire period of the original time series data. Therefore, the process of the above steps corresponds to the above-described division process.
次いでステップS6では、最下位のレイヤーのセグメントに原時系列データの対応部分が取り込まれ、セグメントデータが構成される。 Next, in step S6, the corresponding portion of the original time series data is taken into the segment of the lowest layer, and segment data is configured.
次いでステップS7では、セグメントデータに対して、次式の最大エントロピー法により最小ラグ値から最大ラグ値までバーグの係数Pmを求め、ラグ値に対するPmの推移をステップS8において評価する。
後に例示するように、係数Pmは、ラグ値が増大するに従って単調減少する量であり、ホワイトノイズ的な時系列データでは、この減少は極めてゆるやかであるのに対して、三角関数波に一定量の雑音を重畳した時系列データや、その他の通常観測されるように周期性を有する成分を含む時系列データでは、係数Pmの減少はやや大きくなる。また、周期性を有する成分の他に雑音を含まないか、またはごくわずかに含む三角関数波として記述されるような極めて周期性の高い時系列データでは、係数Pmの減少は極めて急激におこり、最小ラグ値に対する係数Pmの値と、最大ラグ値(セグメントデータのデータ点数−1)に対する係数Pmの比は7〜10桁以上に達する。 As illustrated later, the coefficient P m is an amount that monotonously decreases as the lag value increases. In white noise time-series data, this decrease is very gradual, but is constant in a trigonometric function wave. In time-series data on which a large amount of noise is superimposed and other time-series data including a component having periodicity as normally observed, the decrease in the coefficient P m is slightly increased. In addition, in the time series data with extremely high periodicity described as a trigonometric wave that contains no or very little noise in addition to the periodic component, the coefficient P m decreases very rapidly. The ratio of the coefficient P m to the minimum lag value and the coefficient P m to the maximum lag value (number of data points of segment data−1) reaches 7 to 10 digits or more.
そして、セグメントデータの全期間に渡ってデータ構造が変化せず、均一な時系列データでは、係数Pmは最小のラグ値から増大する際に、最初急激に減少した後に、プラトー領域に至って最大ラグ値まで緩やかに漸減するのに対して、セグメントデータの期間内においてデータ構造が変化する不均一な時系列データでは、プラトー領域において、係数Pmの急激な変化が生じる。 The data structure does not change over the entire period of the segment data. For uniform time-series data, the coefficient P m decreases rapidly from the first when it increases from the minimum lag value, then reaches the plateau region and reaches the maximum. whereas tapers gradually to a lag value, the non-uniform time-series data the data structure is changed in the period of segment data in the plateau region, a rapid change of the coefficients P m occurs.
そこでステップS8においては、セグメントデータの評価として、a:最小ラグ値に対する係数Pmの値と、最大ラグ値に対する係数Pmの比を算出し、その値を設定値と比較して、セグメントデータの周期性の強さを評価することと、b:プラトー領域における急激な値の変化を検出して、データ構造の変化を検出することが行われる。 In Therefore step S8, as the evaluation of the segment data, a: the value of the coefficient P m to the minimum lag value, calculating a ratio of the coefficients P m for the maximum lag value, and compares that value with the set value, the segment data Evaluation of the strength of the periodicity of b, and b: detecting a sudden change in value in the plateau region to detect a change in data structure.
ステップS8においてデータ構造が均一であると判定された場合には、ステップS9に移行して、データ構造を評価したセグメントを基本セグメントとして設定し、次いでステップS10に移行して、全てのセグメントに対してセグメントデータのデータ構造の評価が完了したか否かを判定し、完了していない場合にはステップS6に移行して、次のセグメントに対しての上述した処理が行われる。 If it is determined in step S8 that the data structure is uniform, the process proceeds to step S9, the segment for which the data structure has been evaluated is set as a basic segment, and then the process proceeds to step S10 for all segments. Then, it is determined whether or not the evaluation of the data structure of the segment data is completed. If the evaluation is not completed, the process proceeds to step S6, and the above-described processing is performed for the next segment.
一方、ステップS8においてデータ構造が不均一であると判定された場合には、ステップS11に移行し、ステップS11においては、ステップS8において評価されたセグメントのレイヤーの下位にサブレイヤーを構成し、このサブレイヤーに、データ構造の変化している個所において2分割したセグメントを配置する。 On the other hand, if it is determined in step S8 that the data structure is non-uniform, the process proceeds to step S11. In step S11, a sublayer is configured below the segment layer evaluated in step S8. In the sublayer, a segment divided into two is arranged at a location where the data structure is changed.
次いでステップS6に移行して、ステップS11において構成されたサブレイヤーに配置された2つのセグメントに対して、上述と同様な処理が行われる。即ち、ステップS8において、サブレイヤーのセグメントのデータ構造が均一であると判定された場合には、このセグメントもステップS10において基本セグメントとして設定され、一方、データ構造が不均一であると判定された場合には、サブレイヤーの下位に更にサブレイヤーが構成され、このサブレイヤーに、データ構造の変化している個所において2分割したセグメントが配置され、更にステップS6に移行して上述と同様な処理が行われる。 Next, the process proceeds to step S6, and the same processing as described above is performed on the two segments arranged in the sublayer configured in step S11. That is, if it is determined in step S8 that the data structure of the sub-layer segment is uniform, this segment is also set as a basic segment in step S10, while the data structure is determined to be non-uniform. In this case, a sub-layer is further formed below the sub-layer, and a segment divided into two at the location where the data structure is changed is arranged in this sub-layer, and the process proceeds to step S6 and the same processing as described above. Is done.
こうしてステップS10において、サブレイヤーのセグメントを含む全てのセグメントに対してセグメントデータのデータ構造の評価が完了したか否かを判定し、完了と判定された場合には、ステップS12に移行し、ステップS9において基本セグメントとして設定された全てのセグメントに対応して、記憶手段に記憶領域を設定し、初期化がされて、解析シーケンスが完成する。 Thus, in step S10, it is determined whether or not the evaluation of the data structure of the segment data has been completed for all segments including the sublayer segment. If it is determined that the segment data has been completed, the process proceeds to step S12. Corresponding to all the segments set as basic segments in S9, storage areas are set in the storage means, initialization is performed, and the analysis sequence is completed.
図6は、ある解析対象の時系列データに対して、解析シーケンス決定手段6が決定したレイヤーとセグメントの構成を模式的に示すものであり、レイヤーは、時系列データの全期間データに対応する一つのセグメントが配置されるレイヤー0の下位に、レイヤーN(N:自然数)まで、N+1層のレイヤーが構成されており、また、レイヤーNのセグメント3に対応してサブレイヤー1が構成されている。こうしてレイヤーNの、セグメント3を除く全てのセグメントが基本セグメントとして設定され、またレイヤーNのセグメント3が分割された、サブレイヤー1の2つのセグメントが基本セグメントとして設定されている。そしてこれらの基本セグメントの全てに対応して記憶手段に記憶領域が設定されている。
FIG. 6 schematically shows the structure of layers and segments determined by the analysis sequence determination means 6 for time series data to be analyzed, and the layers correspond to the whole period data of the time series data. N + 1 layers are configured up to layer N (N: natural number) below
データ解析手段7は、解析シーケンス決定手段6により上述したように決定されたレイヤーとセグメントの構成に基づき、セグメントに取り込まれた原時系列データ、即ち、セグメントデータの解析を行う。
The
データ解析手段7によるセグメントデータの解析は、最上位のレイヤー0から最下位のレイヤーN及びサブレイヤー1に向かって順次同様な処理の流れで行われるが、この際、上位のレイヤーのセグメントのセグメントデータは、上述したようにセグメント当たりのデータ点数が予め設定された最大データ点数を越えないように、原時系列データのデータ点数をバンチング処理により低減してセグメントデータを構成して処理が行われる。例えばレイヤー0のセグメントには、バンチング処理により、原時系列データの全期間のデータ点数を1/2Nに低減したデータをセグメントデータとして取り込み、レイヤー1のセグメントには、原時系列データの全期間のデータ点数を1/2N−1に低減したデータをセグメントとして取り込む。
The analysis of the segment data by the data analysis means 7 is performed in the same process flow sequentially from the
更に、上位のレイヤーにおける各セグメントの解析結果中、上位側のセグメントに対して求められたパラメータは、その上位のセグメントに対応する全ての基本セグメントの記憶領域に同様に順次記憶されて行く。一方、上記パラメータから構成される当てはめ曲線と、原時系列データの対応部分との差演算により求められる残差時系列データは次の下位レイヤーのセグメントデータとして設定して解析が行われる。 Further, in the analysis result of each segment in the upper layer, the parameters obtained for the upper segment are sequentially stored in the storage areas of all the basic segments corresponding to the upper segment. On the other hand, the residual time series data obtained by calculating the difference between the fitting curve composed of the above parameters and the corresponding portion of the original time series data is set and analyzed as segment data of the next lower layer.
そこで、データ解析手段7及び解析結果評価手段8における各セグメントの解析処理の流れを図7、図8の流れ図につき説明する。
まずステップS101では、セグメントに原時系列データの対応部分が取り込まれて、セグメントデータが構成される。この際、レイヤー0〜N−1のセグメントの場合には、上述したとおり、バンチング処理によるデータ点数の低減が行われる。
Therefore, the flow of analysis processing of each segment in the data analysis means 7 and the analysis result evaluation means 8 will be described with reference to the flowcharts of FIGS.
First, in step S101, the corresponding portion of the original time series data is taken into the segment to form segment data. At this time, in the case of segments of
次いでステップS102では、ステップS101において取り込まれたセグメントデータに対して、上述のステップS8と同様に最大エントロピー法により最小ラグ値から最大ラグ値までバーグの係数Pmを求め、夫々のラグ値に対するPmの推移を算出する。 Next, at step S102, with respect to the segment data captured in step S101, obtains a coefficient P m of Burg up lag value from the minimum lag value by the maximum entropy method in the same manner as in step S8 described above, P for lag values each Calculate the transition of m .
次いでステップS103では、上述のステップS8と同様に、ラグ値に対するバーグの係数Pmの推移を評価し、それに基づいて最大エントロピー法によるパワースペクトル密度の計算に用いるM個のラグ値を決定する。このラグ値は単数又は複数で良い。 Next, at step S103, similarly to step S8 mentioned above, to evaluate the changes in coefficients P m of Burg for lag value, to determine M lag value used in the calculation of the power spectral density by the maximum entropy method based thereon. This lag value may be singular or plural.
このラグ値の決定は、例えば、a:周期性の強いデータ、即ち最小ラグ値に対する係数Pmの値と最大ラグ値に対する係数Pmの比が大きい場合には、ラグ値を50%以内とし、データ点数が多い場合には、より小さなラグ値とする、b:周期性の弱い場合には、ラグ値は50%以上で、時系列の構造を効率的に取り出すために75%までとし、データ点数が少ない場合には、より大きなラグ値とする、等のラグ値決定のための条件や、その他、数多くの時系列データにおける解析から得られた知見を知識ベースとして記憶しておき、上述したバーグの係数Pmの評価結果から知識ベースを用いて、単数又は複数のラグ値を自動的に決定する構成とすることができる。この他、単一系列解析手段2a,2bのインタラクティブな操作により人がラグ値の決定を行うようにすることができる。更に本発明では、このラグ値決定過程においては、設定される解析の精度に応じて、より高精度の解析において、より大きいラグ値を選択するように知識ベース等により構成されており、高精度の解析は、上述した操作者の操作に加えて、相互相関解析手段3からの再解析指令により設定されるように構成している。 This determination of the lag value is, for example, a: periodic strong data, if i.e. the ratio of the coefficient P m to the value of the coefficient Pm and the maximum lag value for the minimum lag value is large, the lag value is within 50% If the number of data points is large, the lag value should be smaller. B: If the periodicity is weak, the lag value should be 50% or more, and up to 75% in order to extract the time-series structure efficiently. When the number of points is small, the conditions for determining the lag value, such as a larger lag value, and other knowledge obtained from the analysis of many time-series data are stored as a knowledge base, as described above. using a knowledge base from the evaluation results of the coefficient of Berg P m, it can be configured to automatically determine one or more of the lag value. In addition, the person can determine the lag value by interactive operation of the single series analysis means 2a, 2b. Furthermore, in the present invention, the lag value determination process is configured by a knowledge base or the like so as to select a larger lag value in a higher accuracy analysis according to the accuracy of the set analysis. This analysis is configured so as to be set by a reanalysis command from the cross-correlation analysis means 3 in addition to the operation of the operator described above.
なお、最大エントロピー法において時系列が雑音を含まず線スペクトル成分のみからなる特殊な場合を除き、ラグ値をデータ点数−1点まで任意に設定し得て、その際になんら制約の存在しないことは、前掲の非特許文献2において明かにされた重要な知見である。
In the maximum entropy method, except for the special case where the time series does not include noise and consists only of line spectrum components, the lag value can be arbitrarily set up to −1 data points, and there are no restrictions at that time. These are important findings revealed in the aforementioned
次いでステップS103で決定されたM個のラグ値の夫々につき、ステップS104において、上述した文献等に開示される最大エントロピー法によるパワースペクトル密度(MEM-PSD)の計算を行う。 Next, for each of the M lag values determined in step S103, in step S104, the power spectral density (MEM-PSD) is calculated by the maximum entropy method disclosed in the above-described literature.
次いでステップS105では、ステップS104において算出したM個のパワースペクトル密度の夫々から、適数の主要なピークと、そのパワーを抽出する。この主要なピーク及びそのパワーの抽出動作においては、一般に、最初に抽出された主要はピークに近接しているピークを始めとして、より多くのピークを主要ピークとして抽出すると、原時系列データの再現性がより向上する。このため、この抽出過程においては、設定される解析の精度に応じて、より高精度の解析において、より多くのピークを主要ピークとして抽出するように知識ベース等により構成されており、高精度の解析は、操作者の操作に加えて、相互相関解析手段3からの再解析指令により設定されるように構成している。 Next, in step S105, an appropriate number of main peaks and their powers are extracted from each of the M power spectral densities calculated in step S104. In the operation of extracting the main peak and its power, in general, the first extracted main is the peak close to the peak, and more peaks are extracted as the main peak. More improved. For this reason, in this extraction process, it is configured by a knowledge base or the like so that more peaks are extracted as main peaks in higher-precision analysis according to the set accuracy of analysis. The analysis is configured to be set by a reanalysis command from the cross-correlation analysis means 3 in addition to the operation of the operator.
次いでステップS106では、ステップS105で抽出した主要ピークの周波数の逆数を周期の初期値として、非線形最小自乗法により、複数の当てはめ曲線とそのパラメータが算出される。 Next, in step S106, a plurality of fitting curves and their parameters are calculated by the nonlinear least square method using the reciprocal of the main peak frequency extracted in step S105 as the initial value of the period.
次いでステップS107は、ステップS106で求められた当てはめ曲線のパラメータを、初期値と比較して一致度を評価する。この一致度の評価は、例えば、当てはめ曲線の各周期のモードのパワー(振幅の自乗の1/2)を、上記主要ピークのパワーで除した値により評価する。 Next, in step S107, the degree of coincidence is evaluated by comparing the parameters of the fitting curve obtained in step S106 with initial values. For example, the degree of coincidence is evaluated by a value obtained by dividing the power (1/2 of the square of the amplitude) of the mode of each cycle of the fitting curve by the power of the main peak.
このようにしてステップS107において、最も一致度の高い当てはめ曲線を最適当てはめ曲線として選択する。 In this way, in step S107, the fitting curve having the highest degree of coincidence is selected as the optimum fitting curve.
次いでステップS108では、ステップS107において選択された最適当てはめ曲線のパラメータの一致度が、予め設定された一致度以上であるか否かを判定し、設定された一致度以上でない場合には、ステップS111においてラグ値の再設定指示を行ってステップS103に移行し、また予め設定された一致度以上である場合には、ステップS109に移行して、ステップS107において選択された最適当てはめ曲線を、解析結果評価手段8により評価する。 Next, in step S108, it is determined whether or not the degree of coincidence of the parameters of the optimum fitting curve selected in step S107 is equal to or higher than a preset degree of coincidence. In step S103, the process proceeds to step S103. If the degree of coincidence is equal to or higher than the preset matching degree, the process proceeds to step S109, and the optimum fitting curve selected in step S107 is analyzed. Evaluation is performed by the evaluation means 8.
本発明では、この最適当てはめ曲線のパラメータの一致度の判定において、設定される解析の精度に応じて、より高精度の解析において、より高い一致度となるように構成しており、高精度の解析は、相互相関解析手段3からの再解析指令に加えて、操作者の操作によっても設定されるように構成している。 In the present invention, in the determination of the degree of coincidence of the parameters of the optimum fitting curve, the higher degree of coincidence is obtained in the higher accuracy analysis according to the accuracy of the set analysis. The analysis is configured to be set by the operation of the operator in addition to the reanalysis command from the cross correlation analysis means 3.
図7においてステップS109として示す解析結果評価手段8による評価は、図8の流れ図に示す流れで行われる。
まずステップS112では、ステップS107において選出されたパラメータにより、最適当てはめ曲線を再構成し、そのデータ期間に対応する原時系列データとの差の演算を行い、残差曲線を求める。また当てはめ曲線の再現性を評価する諸量、ここでは残差の標準偏差と残差パワーの時間依存が求められる。
The evaluation by the analysis result evaluation means 8 shown as step S109 in FIG. 7 is performed according to the flow shown in the flowchart of FIG.
First, in step S112, an optimal fitting curve is reconstructed using the parameters selected in step S107, and a difference from the original time series data corresponding to the data period is calculated to obtain a residual curve. Also, various values for evaluating the reproducibility of the fitting curve, here, the standard deviation of the residual and the time dependence of the residual power are required.
次いでステップS113においては、ステップS112において求めた当てはめ曲線の再現性を評価する諸量に基づいて、最適当てはめ曲線と、原時系列データの対応部分が十分に一致しているか否かを判定する。 Next, in step S113, it is determined whether or not the optimum fitting curve and the corresponding portion of the original time series data are sufficiently matched based on various amounts for evaluating the reproducibility of the fitting curve obtained in step S112.
ステップS113において十分に一致していると判定された場合には、ステップS114に移行し、求められたパラメータは、そのセグメントに対応する全ての基本セグメントの記憶領域に同様に順次記憶されて行く。一方、上記パラメータから構成される当てはめ曲線と、原時系列データの対応部分との差演算により求められる残差時系列データは次の下位のレイヤーのセグメントデータとして設定された後、ステップS109の処理が終了する。 If it is determined in step S113 that they are sufficiently matched, the process proceeds to step S114, and the obtained parameters are sequentially stored in the storage areas of all the basic segments corresponding to the segment. On the other hand, after the residual time series data obtained by calculating the difference between the fitting curve composed of the above parameters and the corresponding portion of the original time series data is set as the segment data of the next lower layer, the process of step S109 is performed. Ends.
一方、ステップS113において十分に一致していないと判定された場合には、ステップS115に移行し、セグメント解析結果は破棄し、そのセグメントデータを、そのまま次のレイヤーのセグメントに設定した後、ステップS109の処理が終了する。 On the other hand, if it is determined in step S113 that they do not match sufficiently, the process proceeds to step S115, the segment analysis result is discarded, and the segment data is set as it is to the segment of the next layer, and then step S109. This process ends.
このような処理を行うことにより、一括した解析処理では、十分に一致する当てはめ曲線が得られないセグメントデータは、順次、下位の2つのセグメントに分割されて設定され、次は分割された夫々のセグメントデータとして解析処理されるため、一括した解析処理では、データ構造の変動により十分に一致する当てはめ曲線が得られないセグメントデータであっても、最終的にはサブレイヤーの基本セグメントにおいて、十分に一致する当てはめ曲線を求めることができる。 By performing such processing, segment data for which a fitting curve that sufficiently matches cannot be obtained in the batch analysis processing is sequentially divided into two lower segments, and the next segment data is set. Since it is analyzed as segment data, even if it is segment data that does not provide a fitting curve that matches well due to fluctuations in the data structure in the batch analysis process, it will eventually be sufficient in the basic segment of the sublayer. A matching fit curve can be determined.
ステップS109の処理が終了すると、次にステップS110に移行して、全てのセグメントの解析が終了したか否かを判定し、終了していない場合には、ステップS101に移行して、次のセグメントの解析を行い、全てのセグメントの解析処理が終了している場合には、セグメント解析が終了する。 When the process of step S109 is completed, the process proceeds to step S110, where it is determined whether or not the analysis of all segments has been completed. If not, the process proceeds to step S101, and the next segment is determined. If all the segments have been analyzed, the segment analysis is completed.
次に、上述したとおり、一括解析結果構成手段13aは、再構成スペクトル計算手段10とその他計算手段12を利用して図3に示すような解析結果10を構成し、またセグメント解析結果構成手段13bは、再構成スペクトル計算手段10と、頂位位相推移計算手段11と、その他計算手段12を利用して図3に示すような解析結果11を構成するものであり、そこでこれらの計算手段6、7の動作を次に説明する。
Next, as described above, the collective analysis result configuring unit 13a uses the reconstructed
まず図9は一括解析結果構成手段13aにより利用されている場合の再構成スペクトル計算手段10の動作を説明する模式的流れ図であり、まずステップS201では、全ての記憶領域1〜記憶領域2Nに保存されている全てのパラメータ(振幅aj及び周期Tj)と残差時系列データεx(t)を取り出す。尚、取り出されたパラメータの個数は図中に示すように便宜的にL個としている。
First, FIG. 9 is a schematic flowchart for explaining the operation of the reconstructed spectrum calculating means 10 when it is used by the collective analysis result constituting means 13a. First, in step S201, all the
次いで、ステップS202において、時系列データのエネルギーを保存するため、L個のパラメータ、即ち、夫々L個の振幅ajの規格化を行う。この規格化は、解析シーケンス決定手段6よりサブレイヤーが構成されておらず、従って全ての基本セグメントのデータ期間が等しい長さの場合には、全てのモードの振幅ajの自乗を基本セグメントの個数2Nで除して、その平方根を新たな振幅ajとすることで規格化を行うことができる。一方、サブレイヤーが構成されており、従ってデータ期間の長さが異なる基本セグメントが含まれる場合には、全ての基本セグメントについて、対応する記憶域に記された全モードの振幅ajを自乗すると共に、この値に、モードの属する基本セグメントの長さと、全期間の長さの比を乗じ、得られた積の平方根を新たな振幅ajとすることで規格化を行うことができる。 Next, in step S202, in order to preserve the energy of the time-series data, L parameters, that is, L amplitudes a j are normalized. In this normalization, when no sub-layer is configured from the analysis sequence determination means 6, and when the data periods of all the basic segments have the same length, the squares of the amplitudes a j of all the modes are set to the basic segment. Normalization can be performed by dividing the square root by a new amplitude a j by dividing the number by 2N. On the other hand, if sub-layers are configured and therefore basic segments with different data period lengths are included, all basic segments are squared with the amplitudes a j of all modes recorded in the corresponding storage areas. At the same time, normalization can be performed by multiplying this value by the ratio of the length of the basic segment to which the mode belongs and the length of the entire period, and setting the square root of the obtained product as the new amplitude a j .
このようにして振幅が規格化されたL個のパラメータ(振幅aj及び周期Tj)の夫々により、ステップS203において、中心周波数1/Tjにピークをもつ孤立ピーク列Pj、即ち図に示される孤立ピーク列( P1 … Pj … PL )を構成する。
In step S203, the isolated peak train P j having a peak at the
一方、残差時系列データεx(t)については、ステップS204において、最大エントロピー法のスペクトル密度(MEM-PSD)の計算を行い、得られたスペクトル密度と、ステップS203で得られた孤立ピーク列をステップS205において合成することにより、再構成スペクトルが計算される。なお、S204でラグ0とおくと残差時系列のパワーが全周波数帯に均一に分布するスペクトルとなる。 On the other hand, for the residual time series data ε x (t), the spectral density (MEM-PSD) of the maximum entropy method is calculated in step S204, and the isolated spectral peak obtained in step S203 is calculated. By recombining the columns in step S205, a reconstructed spectrum is calculated. If the lag is set to 0 in S204, a spectrum in which the power of the residual time series is uniformly distributed in all frequency bands is obtained.
この際、本発明では、理論形式では有限離散値を扱うコンピュータで構成することができない線スペクトル列となる時系列データのスペクトルを、ピーク周波数から離れるに従って急減する連続関数、例えば余弦函数に微弱な雑音を加えた時系列をラグ2にて最大エントロピー法のスペクトル密度の計算を行って得られたスペクトルをδ函数の原型として利用して再構成することにより、コンピュータを利用して、理論と整合する時系列データの解析を行うことができる。また、これにより、2つの時系列データのクロススペクトルの絶対値と等しくなる、2つの時系列データの積を構成することも可能となる。
At this time, in the present invention, the spectrum of the time series data, which is a line spectrum sequence that cannot be configured by a computer that handles a finite discrete value in the theoretical form, is weak to a continuous function such as a cosine function that rapidly decreases as the distance from the peak frequency increases. Reconstruct the spectrum obtained by calculating the spectral density of the maximum entropy method with
次に図10はセグメント解析結果構成手段13bにより利用されている場合の再構成スペクトル計算手段10と頂位位相推移計算手段11の動作を説明する模式的流れ図である。 Next, FIG. 10 is a schematic flow chart for explaining the operations of the reconstructed spectrum calculating means 10 and the top phase transition calculating means 11 when used by the segment analysis result constituting means 13b.
まずステップS301では、基本セグメントに対応する記憶領域1〜2Nのいずれかが選択され、それに記憶されている、対応する基本セグメントの代表時刻、パラメータ、残差データが取り出される。最初は最も早い代表時刻の記憶領域1の記憶内容が取り出され、次いで後述する繰り返し処理により、順次、時間が経過した基本セグメントに対応する記憶領域の内容が取り出されて、全ての記憶領域の内容が取り出される。
First, in step S301, one of the
取り出された記憶領域のデータ中、パラメータ(振幅aj及び周期Tj)と残差時系列データεx(t)は再構成スペクトル計算手段10に入力され、またパラメータ(周期Tj及び頂位位相φj)は頂位位相推移計算手段11に入力されて処理に供される。 Among the extracted data in the storage area, the parameters (amplitude a j and period T j ) and residual time series data ε x (t) are input to the reconstructed spectrum calculation means 10, and the parameters (period T j and peak position). The phase φ j ) is input to the top phase shift calculation means 11 and used for processing.
再構成スペクトル計算手段10に入力されたパラメータ(振幅aj及び周期Tj)は、一括解析結果構成手段13aにより利用されている場合の上述したステップS202の振幅の規格化の処理は行われず、ステップS203と同様に、ステップS302において、中心周波数1/Tjにピークをもつ孤立ピーク列Pj、即ち図に示される孤立ピーク列( P1 … PL )を構成する。
The parameters (amplitude a j and period T j ) input to the reconstructed
一方、残差時系列データεx(t)については、ステップS303において、最大エントロピー法のスペクトル密度(MEM-PSD)の計算を行い、得られたスペクトル密度と、ステップS302で得られた孤立ピーク列をステップS304において合成することにより、再構成スペクトルが計算される。なお、S303でラグ0とおくと残差時系列のパワーが全周波数帯に均一に分布するスペクトルとなる。 On the other hand, for the residual time series data ε x (t), the spectral density (MEM-PSD) of the maximum entropy method is calculated in step S303, and the obtained spectral density and the isolated peak obtained in step S302. By recombining the columns in step S304, a reconstructed spectrum is calculated. If the lag is set to 0 in S303, a spectrum in which the power of the residual time series is uniformly distributed in all frequency bands is obtained.
一方、頂位位相推移計算手段11に入力されたパラメータ(周期Tj及び頂位位相φj)は、ステップS305において解析されて、周波数と頂位位相との対応関係から、基本セグメントの代表時刻に対する頂位位相の周波数依存が得られる。 On the other hand, the parameters (period T j and top phase φ j ) input to the top phase transition calculating means 11 are analyzed in step S305, and the basic segment representative time is determined from the correspondence between the frequency and the top phase. The frequency dependence of the top phase with respect to is obtained.
これらの得られた解析結果は、ステップS306において所定の記録領域に、基本セグメントの代表時刻の属性として記録される。 These obtained analysis results are recorded in the predetermined recording area as attributes of the representative time of the basic segment in step S306.
次いでステップS307では、全ての基本セグメントに対応する記憶領域のデータが解析されたか否かが判定され、解析されていない記憶領域が残っている場合には、結合子Bを経てステップS301に移行して、次の記憶領域の解析が行われる。一方、ステップS307において、全ての基本セグメントに対応する記憶領域のデータが解析されている場合には、終了処理が行われる。 Next, in step S307, it is determined whether or not the data in the storage area corresponding to all the basic segments has been analyzed. If a storage area that has not been analyzed remains, the process proceeds to step S301 via the connector B. Thus, the next storage area is analyzed. On the other hand, if the data in the storage area corresponding to all the basic segments has been analyzed in step S307, end processing is performed.
次に以上に説明した本発明の解析装置を構成する単一系列解析手段の動作を、具体例につき説明する。
解析する時系列データは、図11に示すもので、データは48時間201点の時系列データである。図から分かるように、24時間モードは、時刻24時においてその位相を反転させている。
Next, the operation of the single series analyzing means constituting the analyzing apparatus of the present invention described above will be described with a specific example.
The time series data to be analyzed is shown in FIG. 11, and the data is time series data of 48 points and 201 points. As can be seen from the figure, in the 24-hour mode, the phase is inverted at 24:00.
まず図11に示される時系列データが解析装置1に入力されると、解析シーケンス決定手段2は上述したように図4、図5に示される処理過程を経て解析シーケンスを決定する。
First, when the time series data shown in FIG. 11 is input to the
この解析装置1では、演算処理の精度上、例えば最大データ点数が300点に設定されているとすると、入力された原時系列データのデータ点数は201点であり、原時系列データのデータ点数が最大データ点数以下であるため、図4のステップS3を経て、レイヤーは一層のみ設定され、このレイヤーに一つのセグメントが配置される。
In the
次いでステップS6において、原系列データの全期間データがセグメントに取り込まれる。次いでステップS7において、最大ラグ値(200点)までバーグの係数Pmの推移が計算される。図12が係数Pmの推移を示すものである。 Next, in step S6, the whole period data of the original series data is taken into the segment. Then, in step S7, the transition of the coefficient P m of Berg are calculated up to the maximum lag value (200 points). Figure 12 illustrates the transition of the coefficient P m.
次いでステップS8において係数Pmの挙動が解析される。まず、最初と最後の係数Pmの比が計算され、その結果、それが10桁を超えるために「周期性の強いデータである」と判断される。 Then the behavior of the coefficients P m is analyzed in step S8. First, the ratio of the first and last coefficients P m is calculated, and as a result, since it exceeds 10 digits, it is determined that the data is “periodic data”.
次いで解析シーケンス決定手段2は、係数Pmの推移を解析し、通常、係数Pmは最初急激に減少したのちに最大ラグ値まで漸減するのに対して、入力された原時系列データでは、ラグ値3〜100(対応するデータ点数は4〜101点)の領域で殆ど係数Pmが減少せず、直後に急激に減少して、その後最大ラグ値ちかくまで再びプラトー領域を形成することから、101点を超える付近でデータ構造が変化していることを検出する。
Then analyzing the
従って解析シーケンス決定手段2は、ステップS11の処理に移行し、サブレイヤーを構成し、データ構造の変化している個所において2分割したセグメントを構成する。 Therefore, the analysis sequence determination means 2 moves to the process of step S11, constitutes a sublayer, and constitutes a segment divided into two at locations where the data structure is changed.
次いでステップS6に移行して、ステップS11において分割されたセグメントの夫々につき、順次ステップS8により、データ構造の均一性を評価する。 Next, the process proceeds to step S6, and the uniformity of the data structure is sequentially evaluated in step S8 for each of the segments divided in step S11.
尚、分割されたサブレイヤーのセグメントの境界にデータが位置しているため、このデータは隣接するセグメントに重複して取り込まれ、結果としてサブレイヤーの各セグメントのデータ点数は、いずれも101点となる。 In addition, since the data is located at the boundary of the segment of the divided sublayer, this data is duplicated and captured in the adjacent segment. As a result, the data score of each segment of the sublayer is 101 points. Become.
まず時間が早い方のセグメントデータについてのステップS8の評価において、係数Pmの推移を計算すると図13に示すように推移し、係数Pmの値はラグ値0から7までで14桁ほど減衰し、上述したようなプラトー領域での急激な減衰領域の存在を検出しない。従って解析シーケンス決定手段2は、ステップS9において、このセグメントを基本セグメントとして設定する。ステップS10からステップS6を経て、時間が遅い方のセグメントデータについても同様な評価がなされ、このセグメントについてもステップS9において基本セグメントとして設定する。
In the evaluation of the step S8 for first time earlier of segment data of the calculation of the transition of the coefficient P m remained as shown in FIG. 13, the value of the coefficient P m is about 14 digits in a
次いで、ステップS10において評価が完了と判定され、ステップS12に移行して、ステップS9において基本セグメントとして設定された全てのセグメントに対応して、記憶手段に記憶領域が設定され、初期化がされて、解析シーケンスが完成する。尚、記憶領域の初期化において、それに対応する基本セグメントがステップS8において「極めて周期性が強い」と判断された場合には、記憶領域にその情報が書き込まれる。 Next, in step S10, it is determined that the evaluation is complete, the process proceeds to step S12, and a storage area is set in the storage means and initialized in correspondence with all the segments set as basic segments in step S9. The analysis sequence is completed. In the initialization of the storage area, when it is determined that the corresponding basic segment is “very periodic” in step S8, the information is written in the storage area.
図14は解析シーケンス決定手段2が決定したレイヤーとセグメントの構成を模式的に示すものであり、レイヤーは、時系列データの全期間に対応する一つのセグメント1が配置されるレイヤー0の下位に、サブレイヤー1が構成され、このサブレイヤー1の2つのセグメントが基本セグメント1、2として設定されている。そしてこれらの基本セグメント1、2に対応して記憶手段に記憶領域1、2が設定されている。
FIG. 14 schematically shows the structure of the layers and segments determined by the analysis sequence determination means 2, and the layers are located below
次いでデータ解析手段3は、図7の流れ図に示す処理過程により各セグメントの解析を行う。この実施例においては、まず最初にレイヤー0のセグメントについてセグメント解析が実行されるが、レイヤー0のセグメント1が含むサブレイヤーの2つの基本セグメント1、2はともに極めて周期性の強いデータであると判別されていることから、このレイヤー0のセグメントの解析は省略することができるが、以降の説明では、レイヤー0のセグメント1についても解析を行うものとしている。
Next, the data analysis means 3 analyzes each segment by the process shown in the flowchart of FIG. In this embodiment, segment analysis is first performed on the segment of
まずデータ解析手段3は、ステップS101においてセグメント1に原時系列データの対応するデータを取り込んでセグメントデータとする。この実施例においては、このセグメント1のデータは、原時系列データそのもので、201点のデータ点数を有している。
First, the data analysis means 3 takes in the data corresponding to the original time series data into the
ここでステップS102において、最大エントロピー法により最小ラグ値から最大ラグ値までバーグの係数Pmを求め、夫々のラグ値に対するPmの推移を算出する。この実施例においては、同様の計算を解析シーケンス決定手段2において行って、図12に示されるような結果を得ているので、本発明の単一系列解析手段2a,2bにおいて、実行された計算を適宜期間保持するように構成することにより、再計算を行う必要がなくなる。 Here in step S102, obtains a coefficient P m of Burg from the minimum lag value to the maximum lag value by the maximum entropy method, we calculate a transition of the P m for lag values each. In this embodiment, the same calculation is performed in the analysis sequence determination means 2 to obtain the result as shown in FIG. 12, so that the calculation executed in the single series analysis means 2a, 2b of the present invention is performed. Is appropriately held for a period, it is not necessary to perform recalculation.
次いでステップS103において、図12に示される推移から特徴的な挙動を示すラグ値を複数選択する。図12において最大の特徴はラグ値2〜101のプラトー領域の存在で、周期性の強いデータなので、データ解析手段3は、自体に構成している知識ベース等を利用した選択手段により、小さめのラグ値を選択することとし、例えばラグ値2・50・101の3つのラグ値を選択し、ステップS104において、これらのラグ値に対応する3つのMEM-PSDを計算する。その結果は、図15に示すとおりであり、左から順にラグ値2・50・101での結果を示す。
Next, in step S103, a plurality of lag values indicating characteristic behavior are selected from the transition shown in FIG. In FIG. 12, the greatest feature is the presence of a plateau region with a lag value of 2 to 101, which is highly periodic data. Therefore, the data analysis means 3 is smaller by a selection means using a knowledge base or the like that is constructed in itself. For example, three lag values of
次いでステップS105では、主要ピークとしてラグ値2では周期の初期値24.66時間(パワー49.80)、ラグ値50では28.06時間(同49.76)、ラグ値101では2モードで34.48時間(同29.81)と18.63時間(同17.31)を抽出し、これを次の非線形最小自乗法による計算の初期値とする。
Next, in step S105, as the main peak, the initial value of the period is 24.66 hours (power 49.80) at
次いでステップS106において、上述した3組の初期値に対して非線形最小自乗法により、夫々の場合の最適当てはめ曲線と曲線パラメータの探索が行われる。この場合、非線形最小自乗法解析においては、初期値の周期の近傍に周期値を探索する領域をあらかじめ設定し、初期値が真の値の近傍にあること、およびその近傍で残差自乗和が唯一の極小値を持つことを期待して計算が行われる。 Next, in step S106, the optimum fitting curve and the curve parameter in each case are searched for by the nonlinear least square method with respect to the above three sets of initial values. In this case, in the non-linear least squares analysis, a region for searching for the period value is set in the vicinity of the period of the initial value in advance, the initial value is in the vicinity of the true value, and the residual sum of squares is in the vicinity. Calculations are performed in the hope of having only one local minimum.
即ち、ステップS106においては、探索過程を監視し、探索中に近傍領域から逸脱する場合や、複数の極小値をもつ場合には、初期値そのものが不適切であると判定して、その結果を使用しない。そこで、ラグ値2と50の初期値は不適切であると判定され、唯一、ラグ値101の初期値に対してのみ計算結果を使用するものとして、最適あてはめ曲線とその曲線パラメータが算出される。図16(a)は算出された最適あてはめ曲線であり、また(b)は原時系列データとの残差曲線を示すものである。 That is, in step S106, the search process is monitored, and when the departure from the neighboring region is made during the search or when there are a plurality of minimum values, it is determined that the initial value itself is inappropriate, and the result is obtained. do not use. Therefore, it is determined that the initial values of the lag values 2 and 50 are inappropriate, and the optimal fitting curve and its curve parameters are calculated only using the calculation result for the initial value of the lag value 101. . FIG. 16A shows the calculated optimum fitting curve, and FIG. 16B shows a residual curve with the original time series data.
次いでステップS107においては、当てはめ曲線のパラメータを、初期値と比較して一致度を評価して、ステップS106で算出した最適あてはめ曲線の妥当性が判定される。即ち、図17は、上記初期値のピーク周波数の逆数とパワーを、最終的に得られた曲線パラメータの周期・振幅とを比較したものである。尚、この場合、図中に示す一致度とは、上述したとおり、モードのパワー(振幅の自乗の1/2)をピークパワーで除した値としている。 Next, in step S107, the matching curve parameter is compared with an initial value to evaluate the degree of coincidence, and the validity of the optimum fitting curve calculated in step S106 is determined. That is, FIG. 17 compares the reciprocal of the peak frequency of the initial value and the power with the period and amplitude of the curve parameter finally obtained. In this case, the degree of coincidence shown in the figure is a value obtained by dividing the mode power (1/2 of the square of the amplitude) by the peak power, as described above.
図17から分かるように、一致度は36時間モードでは0.967、18時間モードでは0.803となる。18時間モードの一致度は余り高くないが、ステップS106において算出された唯一の最適あてはめ曲線とそのパラメータであるので、データ解析手段3は、その結果を、ステップS109、即ち、解析結果評価手段4による評価の処理過程に出力する。 As can be seen from FIG. 17, the degree of coincidence is 0.967 in the 36-hour mode and 0.803 in the 18-hour mode. Although the degree of coincidence in the 18-hour mode is not so high, since it is the only optimum fitting curve calculated in step S106 and its parameter, the data analysis means 3 gives the result to step S109, that is, the analysis result evaluation means 4 Output to the process of evaluation by.
まずステップS112においては、ステップS106〜S108を経て得られた最適あてはめ曲線が再構成され、そのデータ期間に対応する原時系列データとの差の演算を行い、残差曲線を求めると共に、当てはめ曲線の再現性を評価する諸量、ここでは残差の標準偏差と残差パワーの時間依存が求められる。 First, in step S112, the optimum fitting curve obtained through steps S106 to S108 is reconstructed, the difference with the original time series data corresponding to the data period is calculated, a residual curve is obtained, and the fitting curve is obtained. The amount of evaluation of the reproducibility of the error, here, the standard deviation of the residual and the time dependence of the residual power are required.
図18は原時系列データ(点線)と最適あてはめ曲線(実線)とを比較したもの、また図19は残差曲線を、20個の区分領域毎に平均したものである(白丸)。残差の標準偏差は0.836と大きく、また図18から、一見して当てはめが良好でないことが分かる。そして図19から、残差はゼロの回りに一様に分布しておらず、時間依存性が認められる。 FIG. 18 shows a comparison between original time series data (dotted line) and an optimum fitting curve (solid line), and FIG. 19 shows an average of the residual curves for every 20 segmented areas (white circles). The standard deviation of the residual is as large as 0.836, and FIG. 18 shows that the fit is not good at first glance. And from FIG. 19, the residual is not uniformly distributed around zero, and time dependence is recognized.
このような残差の標準偏差の値の大小や、残差の時間依存性等から、データ解析手段3は、ステップS113において、最適当てはめ曲線と、原時系列データの対応部分が十分に一致していないと判定し、ステップS115に移行する。即ち、それまでの解析対象であったレイヤー0のセグメント1のセグメントデータは、そのまま次のレイヤー1のセグメント1、2の解析対象データとして繰り越して設定される。
Due to the magnitude of the standard deviation value of the residual and the time dependency of the residual, the data analysis means 3 sufficiently matches the corresponding portion of the optimal fitting curve and the original time series data in step S113. If it is not determined, the process proceeds to step S115. That is, the segment data of
尚、上述したとおり、原時系列データのデータ点数をバンチング処理により低減してセグメントデータを構成する場合には、下位のレイヤーのセグメントに設定するデータは、必ずしも上位のレイヤーのセグメントデータと同一でなくとも良い。 As described above, when the segment data is configured by reducing the number of data points of the original time series data by bunching processing, the data set for the segment of the lower layer is not necessarily the same as the segment data of the upper layer. Not necessary.
次いでデータ解析手段3は、ステップS109からステップS110を経て、ステップS101に移行し、次にサブレイヤー1の基本セグメント1、2の解析が順次行われる。
Next, the data analysis means 3 goes from step S109 to step S110 to step S101, and then the
まず基本セグメント1については、ステップS101においてセグメントデータとして原時系列データの1〜101点までのデータが取り込まれる。
First, for
次いでステップS102における計算により、図13に示すような係数Pmの推移が得られる。尚、この場合には、数値演算精度のためにラグ値7までの推移が求められる。
Next, the transition of the coefficient Pm as shown in FIG. 13 is obtained by the calculation in step S102. In this case, the transition to the
次いでステップS103では、ラグ値として3と7の二つが決定され、ステップS104において決定されたラグ値でMEM-PSDが計算される。その結果は、図20に示すとおりであり、左から順にラグ値3とラグ値7の結果である。
Next, in step S103, two of 3 and 7 are determined as lag values, and MEM-PSD is calculated with the lag value determined in step S104. The results are as shown in FIG. 20, and are the results of
次いでステップS105では、ラグ値2では初期周期値24.00時間(パワー49.66)を、ラグ値7では初期周期値24.00時間(同50.14)が抽出され、これらの二組の初期値を用いて、ステップS106において非線形最小自乗法により二本の最適あてはめ曲線と曲線パラメータが算出される。この場合には、いずれも周期24.00時間、振幅10.00、頂位位相6.00時が得られる。 Next, in step S105, an initial period value of 24.00 hours (power 49.66) is extracted with a lag value of 2, and an initial period value of 24.00 hours (50.14) is extracted with a lag value of 7, using these two sets of initial values, step S106. The two optimum fitting curves and curve parameters are calculated by the nonlinear least square method. In this case, a period of 24.00 hours, an amplitude of 10.00, and a top phase of 6.00 are obtained.
次いでステップS107では、上述した複数の曲線パラメータの組が比較されるのであるが、この実施例では、いずれのパラメータも一致し、且つ当てはめも良好であるので、その最適当てはめ曲線と、そのパラメータがステップS108を経てステップS109に入力される。 Next, in step S107, the plurality of curve parameter sets described above are compared. In this embodiment, since all the parameters match and the fit is good, the optimum fit curve and its parameters are It inputs into step S109 through step S108.
ステップS112において得られる残差の標準偏差は、10−10と十分に小さく、且つ、残差のパワー推移に時間依存がないため、ステップS113において、今度は、最適当てはめ曲線と、原時系列データの対応部分が十分に一致していると判定し、ステップS114に移行する。 Since the standard deviation of the residual obtained in step S112 is sufficiently small as 10−10 and the power transition of the residual does not depend on time, in step S113, the optimum fitting curve and the original time series data are now used. Are determined to match sufficiently, and the process proceeds to step S114.
ステップS114においては、解析結果のパラメータを記憶領域1に保存する。尚、この実施例と異なり、下位に更にサブレイヤーがある場合には、残差が計算されて、それを下位のサブレイヤーに設定するべく保存がされる。
In step S114, the parameter of the analysis result is stored in the
次いでステップS110を経て、ステップS101に移行して、次の基本セグメント2についても同様な解析処理が行われ、こうして、記憶領域1、2に、原時系列データを一般化三角多項式で記述し、且つ、以降の解析に必要な諸量の全てが蓄積される。
Subsequently, the process proceeds to step S101 through step S110, and the same analysis processing is performed for the next
そして記憶領域1、2に蓄積されている原時系列データの解析に必要な諸量は、上述したとおり、一括解析結果構成手段13aまたはセグメント解析結果構成手段13bにより処理されて、図3に示すような解析結果14と15を構成することができる。
As described above, various quantities necessary for analyzing the original time series data stored in the
一例として、一括解析結果構成手段13aにより利用されている場合の再構成スペクトル計算手段6による再構成スペクトルの算出の流れを説明する。
まず図9において、ステップS201では、記憶領域1と2からパラメータ(一般化三角多項式のパラメータ)を取り出される。この際、取り出されたパラメータは2つで、いずれも周期24時・振幅10・頂位位相6時である。
As an example, the flow of calculation of a reconstructed spectrum by the reconstructed spectrum calculating means 6 when used by the collective analysis result configuring means 13a will be described.
First, in FIG. 9, in step S <b> 201, parameters (generalized triangular polynomial parameters) are extracted from the
次いでステップS202で振幅が規格化される。この規格化は、2つのパラメータそれぞれについて、振幅の自乗=100を基本セグメント数(=2)で割って、その平方根(=√50)をあらたな振幅とすることで行われる。 Next, in step S202, the amplitude is normalized. This normalization is performed by dividing the square of amplitude = 100 by the number of basic segments (= 2) and setting the square root (= √50) to a new amplitude for each of the two parameters.
次いでステップS203における孤立ピークのPSDの計算により、同じピーク周波数(=1/24時間)をもつ二つのピークスペクトルが得られる。この孤立ピークの計算は、ピーク周波数の逆数の周期をもつ三角関数のMEM-PSDを求めておき、その面積を各モードのパワーに一致させる規格化により実現する。図21は余弦関数のMEM-PSDの例である。 Next, two peak spectra having the same peak frequency (= 1/24 hours) are obtained by calculating the PSD of the isolated peak in step S203. The calculation of the isolated peak is realized by obtaining a trigonometric MEM-PSD having a period of the reciprocal of the peak frequency and normalizing the area to match the power of each mode. FIG. 21 shows an example of a cosine function MEM-PSD.
このように本発明では、孤立ピーク列は、δ函数の原型を用いて全周波数帯に裾の広がる夫々一つのスペクトルピークを構成することにより、パラメータから得られる各周波数に対応する孤立ピーク群を構成するもので、こうして理論形式では有限離散値を扱うコンピュータで構成することができない線スペクトル列となる時系列データのスペクトルを、ピーク周波数から離れるに従って急減する連続関数、例えば余弦函数に微弱な雑音を加えた時系列をラグ2にて最大エントロピー法のスペクトル密度の計算を行って得られたスペクトルをδ函数の原型として利用して再構成することにより、コンピュータを利用して、理論と整合する原時系列データの再構成スペクトルを求めることができる。尚、この手法は、後述するように、相互相関解析手段3における計算においても同様に用いられる。
As described above, in the present invention, the isolated peak sequence is formed by using one prototype of the δ function to form one spectral peak having a base extending over the entire frequency band, thereby obtaining an isolated peak group corresponding to each frequency obtained from the parameters. In this way, the spectrum of time-series data, which is a line spectrum sequence that cannot be configured by a computer that handles finite discrete values in the theoretical form, is a continuous function that suddenly decreases as it goes away from the peak frequency, for example, noise that is weak to cosine function By using the computer to reconstruct the spectrum obtained by calculating the spectral density of the maximum entropy method at
一方、各記憶領域1、2から取り出された基本セグメント1、2の残差成分から全測定長にわたる残差曲線を求め、次いでステップS204において、残差のMEM-PSDが計算される。次いで、ステップS205において、孤立ピーク群のスペクトルと残差のMEM-PSDとが重畳されて、最終的に原時系列データの一般化三角多項式表現と整合する再構成スペクトルが得られる。なお、S204でラグ0とおくと残差時系列のパワーが全周波数帯に均一に分布するスペクトルとなる。
On the other hand, a residual curve over the entire measurement length is obtained from the residual components of the
相互相関解析手段3の実施の形態を説明すると、図22に示す相互相関解析手段3は、一般化三角多項式パラメータ16a,16bからパワースペクトル密度を求めるパワースペクトル密度計算手段18a,18bと、クロススペクトル(コスペクトル、クオドスペクトル)とコヒーレンスを含む相互相関を定量化する諸量20を求める相互相関諸量計算手段17と、相互相関諸量計算手段17により求められたコヒーレンスを理論値と比較して、それらの差が所定値以下の場合に、その差をなくすようにクロススペクトルを補正する整合性評価補正手段21と、整合性評価補正手段21において評価した差が所定値以上の場合に、単一系列解析手段2a,2bに再解析の指令23a,23bを出力し、所定値以下の場合には、相互相関諸量計算手段17により計算された諸量20と、パワースペクトル密度計算手段18a,18bにより計算されたパワースペクトル密度19a,19b及び整合性評価補正手段により補正されたクロススペクトルを出力する分岐出力手段22を構成したものである。
An embodiment of the cross-correlation analysis means 3 will be described. The cross-correlation analysis means 3 shown in FIG. 22 includes power spectrum density calculation means 18a and 18b for obtaining power spectrum density from the generalized triangular
まず、相互相関解析手段3に入力される2系列の時系列データは、上述と同様に夫々、次式の一般化三角多項式で表現される。なお、クロススペクトルに関わる次式以下の展開は後述の非特許文献3で詳述される。
上式のx(t)とy(t)は、上述したとおり有限長の時系列データを記述するものであるが、仮に、その記述を±∞に延長したとき、それら延長時系列データのパワースペクトル密度Sxx(f)とSyy(f)は次式で表される。
即ち、パワースペクトル密度計算手段18a,18bにおいては、上式に基づき、夫々の時系列データのパワースペクトル密度Sxx(f)とSyy(f)が求まる。この際には、上述したように、線スペクトル列となる時系列データのスペクトルを、ピーク周波数から離れるに従って急減する連続関数、例えば余弦函数に微弱な雑音を加えた時系列をラグ2にて最大エントロピー法のスペクトル密度の計算を行って得られたスペクトルをδ函数の原型として利用して再構成することにより、有限離散値を扱うコンピュータを利用して計算を行うことができる。このδ函数の原型は全周波数帯で正の値をもつ。
That is, in the power spectrum density calculation means 18a and 18b, the power spectrum densities S xx (f) and S yy (f) of each time series data are obtained based on the above equation. In this case, as described above, the time series data spectrum that becomes the line spectrum sequence is continuously reduced at the
すなわち上式はδ函数の原型を用いて次式となる。
一方、相互相関を定量化する諸量のうち、そのクロススペクトルやコヒーレンス等は次式で表される。
同様に回転スペクトルについても、例えば反時計回りの回転スペクトルS+(f)と回転係数CR(f)は次式で表される。
以上の相互相関を定量化する諸量20は、相互相関諸量計算手段17において計算される。
The
本発明では、相互相関解析手段3の整合性評価補正手段21において、計算結果と理論上期待される値との一致を評価する。具体的には、まず、パワースペクトル密度計算手段18a,18bにおいて求められたパワースペクトル密度Sxx(f),Syy(f)と、相互相関諸量計算手段17において求められたクロススペクトルKxy(f),Qxy(f)とから求められるコヒーレンスcoh2(f)を理論値1と比較し、それらの差が所定値以下の場合に整合性が高いと判定する。
In the present invention, the consistency evaluation correction means 21 of the cross correlation analysis means 3 evaluates the coincidence between the calculation result and the theoretically expected value. Specifically, first, the power spectrum density S xx (f), S yy (f) obtained by the power spectrum density calculating means 18 a, 18 b and the cross spectrum K xy obtained by the cross correlation
なお、理論上コヒーレンスは全周波数帯で1となる。この重要な知見は、2008年12月25日刊行予定の下記非特許文献3において証明が記載される。すなわち、コヒーレンスは0〜1の間で連続的に変動するという従来の理解が誤りであったことが明らかとされる。
判定の結果、整合性が高くないと判定された場合には、それらの差を判定し、その差が微細と判定される場合は、整合性評価補正手段21において、その差をなくすようにクロススペクトルを補正する。具体的には、コヒーレンスcoh2(f),フェイズθxy(f),時間遅れτxy(f)、パワースペクトル密度Sxx(f),Syy(f)を固定して、物理量の間の定義式を満足するようにクロススペクトルKxy(f),Qxy(f)を補正する。 As a result of the determination, if it is determined that the consistency is not high, the difference between them is determined. If the difference is determined to be fine, the consistency evaluation correction means 21 crosses the difference so as to eliminate the difference. Correct the spectrum. Specifically, coherence coh 2 (f), phase θ xy (f), time delay τ xy (f), power spectral density S xx (f), S yy (f) are fixed, and between physical quantities The cross spectrums K xy (f) and Q xy (f) are corrected so as to satisfy the definition formula.
このようにして判定の結果、整合性が高いと判定された場合や、差が微細で補正がされた場合には、分岐出力手段22により、相互相関諸量計算手段17により計算された諸量と、パワースペクトル密度計算手段18a,18bにより計算されたパワースペクトル密度Sxx(f),Syy(f)及び整合性評価補正手段21により補正されたコスペクトルKxy(f)、クオドスペクトルQxy(f)が出力される。
As a result of the determination, when it is determined that the consistency is high, or when the difference is fine and corrected, the quantities calculated by the cross-correlation quantity calculation means 17 by the branch output means 22. And the power spectrum density S xx (f), S yy (f) calculated by the power spectrum density calculating means 18a, 18b, the cospectrum K xy (f) corrected by the consistency
一方、整合性が高くないと判定された場合には、以上の計算結果は破棄されると共に、分岐出力手段22を介して、夫々の単一系列解析手段2a,2bに上述した再解析指令を発する。 On the other hand, if it is determined that the consistency is not high, the above calculation result is discarded, and the reanalysis command described above is sent to each single series analysis means 2a, 2b via the branch output means 22. To emit.
こうして本発明では、相互相関解析手段3において、一般化三角多項式パラメータにより求められたコヒーレンスを理論値と比較して整合性を評価し、所定以上の整合性が得られない場合には、単一系列解析手段2a,2bに対して再解析指令を発して、より高精度の解析を行わせるという、フィードバックと繰り返し動作を有するため、2系列の時系列データの間の相互相関を定量化する諸量を非常に高精度に求めて出力することができる。
Thus, in the present invention, the
以上説明したように、本発明においては、2系列の時系列データの夫々を単一系列解析手段2a,2bにおいて、まず最大エントロピー法を適用して高精度なパワースペクトル密度を求めると共に、その顕著なピークから、一般化三角多項式を解く初期値としての、項数と、各三角項の周期値を求め、次いで、これらの初期値により非線形最小自乗法を用いて、各三角項のパラメータを求めると共に、残差時系列データを求め、これらを加算して時系列データを再構成するという解析を行うことにより、何らの数値モデルも想定せず、夫々の時系列データのみを用いて高精度の解析を行え、これらの高精度の解析によって得られた2系列の時系列データの一般化三角多項式パラメータにより、相互相関解析手段3において、2系列の時系列データの間の相互相関を定量化する諸量を求めることができる。 As described above, in the present invention, the single-sequence analysis means 2a, 2b first applies a maximum entropy method to obtain a highly accurate power spectral density for each of two series of time series data, From the peak, find the number of terms as the initial value for solving the generalized trigonometric polynomial and the periodic value of each triangular term, and then obtain the parameters of each triangular term using these initial values using the nonlinear least squares method In addition, by calculating the residual time series data and adding them to reconstruct the time series data, it is possible to obtain a high-precision data using only each time series data without assuming any numerical model. The cross-correlation analysis means 3 uses the generalized trigonometric polynomial parameters of the two series of time series data obtained by these high precision analyses. It can be determined quantities to quantify the cross-correlation between the data.
また本発明においては、2系列の時系列データの間の相互相関を定量化する諸量を求める他、予測(および解析)のための数値モデルを全く要さず、2系列の時系列データのみを用いて予測支援機能を発揮させることができる。 In the present invention, in addition to obtaining various quantities for quantifying the cross-correlation between two series of time series data, no numerical model for prediction (and analysis) is required, and only two series of time series data are used. Can be used to demonstrate the prediction support function.
即ち、図23は本発明の装置を利用して構成した予測支援機能を模式的に示すもので、実線の矢印はデータの流れを示し、破線の矢印は制御(フィードバック)の流れを示すものである。 That is, FIG. 23 schematically shows a prediction support function configured using the apparatus of the present invention, where solid arrows indicate the flow of data and broken arrows indicate the flow of control (feedback). is there.
なお、予測支援機能の目的は、多数の時系列データを解析してその中に、ある時系列Aの将来挙動が別の時系列Bの現在までの挙動により予想し得る一対の時系列AとBを見いだし、現在のBによりAの将来を予測することである。 The purpose of the prediction support function is to analyze a large number of time series data and include a pair of time series A in which the future behavior of one time series A can be predicted by the behavior of another time series B to the present. Find B and predict the future of A with the current B.
まず、予測支援機能の動作の概要を説明すると、図中符号24a,24b,24c,24d,…は複数の系列1,2,3,4,…の時系列データと、本発明の装置による解析結果の諸量を示すもので、予測対象系列の参照系列の決定処理手段25においては、多数の時系列とその解析結果を検査して共通の周期をもつ有用なモード(指標モード)をともに含む一対の時系列(予測対象系列と参照系列)を抽出し、次いで算出処理手段32において、基底変動や指標モードの解析区間と予測区間における挙動など、予測支援に直接用いる諸量を求めて、それらの諸量33を出力する。
First, the outline of the operation of the prediction support function will be described. In the figure,
上記動作を更に詳細に説明すると、まず処理手段26においては、同じ周期の三角函数波を含む一対の時系列データを探索する(それらを以下指標モードとする)。指標モードとしては、(a)長周期であること、(b)指標モードの振幅がその他のモードのそれに対して十分な大きさのあること、(c)頂位位相が接近または一方に半周期加えた位相が接近していること(これは必須ではない)、等が求められる。こうして位相が先行する指標を参照指標モード28、後行する指標を予測対象指標モード27として求める。
The above operation will be described in more detail. First, the processing means 26 searches for a pair of time series data including triangular function waves having the same period (hereinafter referred to as index mode). The index mode includes (a) a long period, (b) the amplitude of the index mode is sufficiently large compared to that of the other modes, and (c) the top phase approaches or half of the period It is required that the added phase is close (this is not essential). In this way, the index with the leading phase is obtained as the
そして本発明では、参照指標モードの現在の挙動から予測対象指標モードの将来の挙動を推測する。このため、それぞれの時系列について指標モードが時系列全体の挙動とその他の各モードの挙動に対して重要な立ち位置を占めていることが肝要である。処理手段29においては、処理手段26において見いだした指標モードをこの視点で評価し、図中右側の破線矢印に示すように、不適切ならば同じ時系列の組の別の指標モードを探索するか、別の時系列の組の探索に移る。 In the present invention, the future behavior of the prediction target index mode is estimated from the current behavior of the reference index mode. For this reason, it is important that the index mode occupies an important position for the behavior of the entire time series and the behavior of each other mode for each time series. In the processing means 29, the index mode found in the processing means 26 is evaluated from this viewpoint, and as shown by the broken line arrow on the right side in the figure, if it is inappropriate, another index mode of the same time series set is searched. Move on to search for another time series set.
更に処理手段29においては、指標モードよりも長周期のモードの一群と短周期のモードの一群について、指標モードとの兼ね合いが検討される。(a)長周期モードの一群については、極長周期の少数個のモードで時系列全体の基本的挙動(基底変動)を記述できることが望ましい。なぜならば、そのような長周期モードは解析区間の外側、予測区間においても大きく変動しないことが期待されるからである。また、(b)短周期モードの一群については、時系列の挙動への寄与が極力少ないことが望ましい。例えば一週間のリズムなどの外部要因によることが明らかなモードを除いては、それらの寄与が大きい場合、予測区間において予測対象指標モードが期待通りに推移しても、全モードの和・実際の観測値の推移との間に大きな乖離が予想されるからである。 Further, the processing means 29 examines the balance between the index mode for the group of long-period modes and the group of short-period modes than the index mode. (a) For a group of long-period modes, it is desirable that the basic behavior (basic variation) of the entire time series can be described with a small number of extremely long-period modes. This is because such a long period mode is expected not to fluctuate greatly in the prediction section outside the analysis section. In addition, it is desirable that (b) a group of short period modes contributes as little as possible to the time-series behavior. Except for modes that are apparently due to external factors such as the rhythm of a week, if the contribution is large, even if the prediction target indicator mode changes as expected in the prediction interval, the sum of all modes and the actual This is because a large divergence from the transition of the observed values is expected.
以上の処理により、予測支援に用いる一対の時系列と解析結果30、31として得ることができ、次いで処理手段32においては、こうして得られた2つの時系列から、実際に予測支援に用いるさまざまな値、即ち上述した諸量33を得るための処理が行われる。この処理においても、結果が不満足な場合にはその一対の時系列は予測支援に適さないとして、図中左側の破線矢印で示すように処理手段25に指令がフィードバックされる。
Through the above processing, a pair of time series used for prediction support and
処理手段32においては、(1)に示すように、それぞれの時系列について、その三角函数波のうちで基底変動を構成するモードが抽出され、それに指標モードが加えられる。また(2)に示すように、少数個の長周期モードで観測期間全体にわたる時系列の基本的挙動を再現できるかが調べられ、できない場合には、時系列(の対)は破棄され、処理手段25における処理に移行する。 In the processing means 32, as shown in (1), for each time series, a mode constituting a base variation is extracted from the triangular function wave, and an index mode is added thereto. In addition, as shown in (2), it is checked whether the basic behavior of the time series over the entire observation period can be reproduced in a small number of long-period modes. If not, the time series (pair) is discarded and processed. The process proceeds to means 25.
また、処理手段32においては、(3)に示すように、二つの指標モードがそれぞれの時系列の解析区間全体において安定しているかが調べられる。即ち、解析区間の前半・真ん中・後半といった具合にデータを複数のセグメントに分割し、その各セグメントにおいて、与えられた周期値をもつ一般化三角多項式でデータに当てはめを行い、いずれのセグメントでも少なくとも指標モードが安定して存在するか、または変化がゆるやかであることを確認する。そして、加えて、(4)に示すよう、2つの指標モードがそれぞれの時系列の解析区間全体において安定した寄与を示すことを確認する。これは例えば、指標モードの一周期ごとに時系列の平均値を求め、その推移がゆるやかであるかなどを確認しておこなう。更に、その他、時系列に応じて(外部からの強制周期の有無、利用可能な付加情報の有無など)、二つの系列を予測対象系列・参照系列として扱って良いかが検討される。 Further, in the processing means 32, as shown in (3), it is checked whether the two index modes are stable in the entire analysis section of each time series. That is, the data is divided into a plurality of segments such as the first half, the middle, and the second half of the analysis interval, and each segment is fitted to the data with a generalized trigonometric polynomial having a given period value. Make sure that the indicator mode exists stably or changes slowly. In addition, as shown in (4), it is confirmed that the two index modes show a stable contribution in the entire analysis section of each time series. This is performed, for example, by obtaining an average value of time series for each cycle of the index mode and confirming whether the transition is slow. In addition, depending on the time series (existence of forced period from outside, presence / absence of available additional information, etc.), it is examined whether the two series can be handled as the prediction target series / reference series.
次に本発明に係る2系列の時系列データの解析装置による解析例を説明する。 Next, an analysis example by the two-series time-series data analysis apparatus according to the present invention will be described.
まず図24は、共に6つの正弦波に雑音を重畳した次式で示される2つの時系列を解析した例である。
そして図中(c)はクロススペクトルの絶対値、(d)はコスペクトル、(e)は頂位時刻の差(時間遅れ)、および(f)はクオドスペクトルの解析結果である。なお、時間遅れは余弦函数の頂位時刻について計算されている。
上述したように与えたパラメータから、共通の周期である2、5、7、13、すなわち周波数1/2、1/5、1/7、1/13で相互相関を持つことが期待されるが、図の(c)〜(f)ではまさにこの4つの周波数近傍に相互相関が認められる。
図中の(g)に、解析して得たそれぞれの最適あてはめ曲線のパラメータと時間遅れを示しているが、パラメータは時系列を発生する際に設定した値とよく一致していることが分かる。また時間遅れについても、設定値から期待される値によく一致している。大きな雑音の重畳にもかかわらず理論値と一致する結果が得られたことは、MEM-PSDのピークとピークパワーを初期値とする非線形最小自乗法を解くという、本発明に係る解析の手法が正しく機能しいてる証左である。
First, FIG. 24 shows an example in which two time series represented by the following equation in which noise is superimposed on six sine waves are analyzed.
In the figure, (c) is the absolute value of the cross spectrum, (d) is the cospectrum, (e) is the difference in peak time (time delay), and (f) is the analysis result of the quad spectrum. The time delay is calculated for the crest function's top time.
From the parameters given as described above, it is expected that there is a cross-correlation at 2, 5, 7, 13 that is a common period, that is,
(G) in the figure shows the parameters and time delay of each optimum fitting curve obtained by analysis, but it can be seen that the parameters agree well with the values set when generating the time series. . The time delay also agrees well with the value expected from the set value. The fact that the result agrees with the theoretical value despite the large noise superposition is that the analysis method according to the present invention, which solves the nonlinear least square method with the peak and peak power of the MEM-PSD as the initial value, is used. It is proof that it is functioning correctly.
次に図25は、図中(a)に示す振幅変調する時系列データと、(b)に示す正弦波時系列データの2つの時系列を解析した例であり、夫々の時系列データは次式により示される。
上式より、x(t)は周波数0.5に主たるピークを、周波数0.25と0.75に、振幅変調に由来するピークを持つ。また、y(t)は周波数0.5にピークを持つ。そして相互相関は周波数0.5においてのみ期待される。
図24と同様に(c)はクロススペクトルの絶対値、(d)はコスペクトル、(e)は頂位時刻の差(時間遅れ)、および(f)はクオドスペクトルの解析結果である。
この解析で得たクロススペクトルの絶対値(=二つのスペクトルの積の平方根)では、振幅変調に由来する周波数0.25と0.75の位置を含め、3つのピークが認められる。これは解析において、孤立ピークを全周波数帯にわたり構成した結果である(全周波数帯に対して構成しなければ振幅変調に由来するピークはつぶれてしまう)。
また、その他の3つの図、(d)のコスペクトル、(f)のクオドスペクトル及び(e)の頂位時刻の差(=時間遅れ)では、いずれも周波数0.5においてのみ、相互相関が認められる。
図中(g)に示すように、解析で得られたそれぞれの最適あてはめ曲線のパラメータと時間遅れの値は、いずれも設定したパラメータ値をよく再現していることが分かる。
Next, FIG. 25 is an example in which two time series of time series data for amplitude modulation shown in (a) and sine wave time series data shown in (b) are analyzed. It is shown by the formula.
From the above equation, x (t) has a main peak at frequency 0.5, and peaks at frequencies 0.25 and 0.75 due to amplitude modulation. Y (t) has a peak at a frequency of 0.5. And cross-correlation is only expected at frequency 0.5.
As in FIG. 24, (c) is the absolute value of the cross spectrum, (d) is the co-spectrum, (e) is the difference in peak time (time delay), and (f) is the analysis result of the quad spectrum.
In the absolute value of the cross spectrum obtained by this analysis (= square root of the product of the two spectra), three peaks including the positions of the frequencies 0.25 and 0.75 derived from the amplitude modulation are recognized. This is a result of configuring isolated peaks over the entire frequency band in the analysis (if not configured for the entire frequency band, the peak derived from amplitude modulation is crushed).
In the other three figures, (d) co-spectrum, (f) quad spectrum, and (e) peak time difference (= time lag), cross-correlation is recognized only at frequency 0.5. .
As shown in (g) in the figure, it can be seen that the parameters of the optimum fitting curve and the time delay values obtained by the analysis both closely reproduce the set parameter values.
次に図26は、次式で示される図中(a)のうなりの時系列x(t)と、周期的方形波y(t)の2つの時系列を解析した例である。
計算式より、周波数1/50においてのみ、相互相関が期待される。解析結果を示す(c)のクロススペクトルの絶対値はこの位置にピークを持つと共に、周期的方形波の高周波モードの位置にピークを持つ。ここで、注目すべきは、うなりのスペクトルは周波数1/50と1/51の近接ピークをもつが(非特許文献2、pp.198-199、または解析例3の図下段の数表参照)、クロススペクトルの絶対値の計算結果には理論式から期待される周波数1/50のピークのみが反映されていることである。本解析システムの周波数分解能の高さの証左である。
図中(g)に示すように、解析で得られたそれぞれの最適あてはめ曲線のパラメータと時間遅れの値は、いずれも設定したパラメータ値をよく再現していることが分かる。ここに時間遅れは+1.135となり、理論値0とは、周期に対する誤差約2%で一致している。
Next, FIG. 26 is an example in which two time series of a beat time series x (t) of (a) and a periodic square wave y (t) in the diagram expressed by the following equation are analyzed.
From the calculation formula, cross-correlation is expected only at
As shown in (g) in the figure, it can be seen that the parameters of the optimum fitting curve and the time delay values obtained by the analysis both closely reproduce the set parameter values. Here, the time delay is +1.135, which is the same as the
以上に説明したように、本発明においては、2系列の時系列データの夫々を単一系列解析手段において、最大エントロピー法と非線形最小自乗法を用いて、各三角項のパラメータを求めると共に、残差時系列データを求め、これらを加算して時系列データを再構成するという解析を行うことにより、何らの数値モデルも想定せず、夫々の時系列データのみを用いて高精度の解析を行え、これらの高精度の解析によって得られた2系列の時系列データの一般化三角多項式パラメータにより、相互相関解析手段において、2系列の時系列データの間の相互相関を定量化する諸量を求めることができる As described above, in the present invention, the parameters of each triangular term are obtained by using the maximum entropy method and the nonlinear least square method in the single sequence analysis means for each of the two series of time series data. By calculating the difference time series data and adding them to reconstruct the time series data, it is possible to perform high-precision analysis using only each time series data without assuming any numerical model. Based on the generalized trigonometric polynomial parameters of the two series of time series data obtained by these high-accuracy analyses, various quantities for quantifying the cross correlation between the two series of time series data are obtained in the cross correlation analysis means. be able to
本発明に係る解析装置においては、実在の時系列データに対して、時間軸上での一般化三角多項式表現と周波数軸上での再構成スペクトル表現を機械的に計算する。このとき、時系列データ・一般化三角多項式表現・再構成スペクトル表現それぞれの間には、コンピュータ上に構築された数値演算処理システムの能力の制限の元で、理論上の整合性が完全に保証される。すなわち、任意の実データに対して時間軸上および周波数軸上のもっとも整合的で透徹した理解が機械的に得られる。 The analysis apparatus according to the present invention mechanically calculates a generalized trigonometric polynomial expression on the time axis and a reconstructed spectrum expression on the frequency axis for real time series data. At this time, the theoretical consistency between the time series data, generalized trigonometric polynomial expression, and reconstructed spectrum expression is completely guaranteed under the limitations of the numerical processing system built on the computer. Is done. That is, the most consistent and transparent understanding on the time axis and the frequency axis can be mechanically obtained for arbitrary actual data.
また、本発明に係る解析装置では、何らの数値モデルも想定せず、時系列データのみを用いて高精度の解析を行え、しかもその解析は、通常のパーソナルコンピュータ程度の演算能力でも十分に実用的であるという効果があり、特に、解析対象の時系列データを、そのデータ点数やデータ構造に応じて適切にセグメント化して処理を行うことにより、解析する時系列データによらず、常に適切に高精度な解析を行えるという格別なる効果がある。 In addition, the analysis apparatus according to the present invention can perform high-precision analysis using only time-series data without assuming any numerical model, and the analysis is sufficiently practical even with a computation capability comparable to that of a normal personal computer. In particular, the time series data to be analyzed is appropriately segmented according to the number of data points and the data structure and processed, so that it is always appropriate regardless of the time series data to be analyzed. There is a special effect that high-precision analysis can be performed.
そして、本発明は、実在の時系列データ、即ち、ある限られた期間に観測されたある点数の有限長離散時系列のデータに加えて、例えば地層断面の距離に対する輝度変化等の、空間的に連続するデータの解析にも利用できるものである。 In addition to real time-series data, that is, finite-length discrete time-series data of a certain number of points observed in a certain limited period, the present invention provides spatial data such as a change in brightness with respect to the distance of the geological section. It can also be used for analysis of continuous data.
本発明を利用して構成される予測支援システムは、多数の時系列中から一対の時系列を見いだすことにより、一方の時系列の現在までの挙動から他方の時系列の将来挙動を予測することが期待される。 The prediction support system configured by using the present invention predicts the future behavior of one time series from the current behavior of one time series by finding a pair of time series from among many time series. There is expected.
1 統合構造解析手段
2a,2b 単一系列解析手段
3 相互相関解析手段
4a,4b 時系列データ
5a,5b,5c 諸量
6 解析シーケンス決定手段
7 データ解析手段
8 解析結果評価手段
9 制御手段
10 再構成スペクトル計算手段
11 頂位位相推移計算手段
12 その他計算手段
13a 一括解析結果構成手段
13b セグメント解析結果構成手段
14 一括解析結果
15 セグメント解析結果
16a,16b 一般化三角多項式パラメータ
17 相互相関諸量計算手段
18a,18b パワースペクトル密度計算手段
19a,19b パワースペクトル密度
20 諸量
21 整合性評価補正手段
22 分岐出力手段
23a,23b 再解析の指令
24a,24b,24c,24d 諸量
25 決定処理手段
26 処理手段
27 予測対象指標モード
28 参照指標モード
29 処理手段
30,31 解析結果
32 処理手段
33 諸量
DESCRIPTION OF
Claims (26)
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2008282192A JP4711356B2 (en) | 2008-10-31 | 2008-10-31 | Two-series time-series data analysis apparatus and computer-readable recording medium having recorded two-series time-series data analysis program |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
JP2008282192A JP4711356B2 (en) | 2008-10-31 | 2008-10-31 | Two-series time-series data analysis apparatus and computer-readable recording medium having recorded two-series time-series data analysis program |
Publications (2)
Publication Number | Publication Date |
---|---|
JP2010108402A true JP2010108402A (en) | 2010-05-13 |
JP4711356B2 JP4711356B2 (en) | 2011-06-29 |
Family
ID=42297761
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
JP2008282192A Expired - Fee Related JP4711356B2 (en) | 2008-10-31 | 2008-10-31 | Two-series time-series data analysis apparatus and computer-readable recording medium having recorded two-series time-series data analysis program |
Country Status (1)
Country | Link |
---|---|
JP (1) | JP4711356B2 (en) |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107219840A (en) * | 2017-05-05 | 2017-09-29 | 浙江理工大学 | Towards the regulating valve nonlinear characteristic detection method and system of Natural Gas Station |
JP2019191630A (en) * | 2018-04-18 | 2019-10-31 | 国立大学法人京都大学 | Analysis device, analysis method and computer program |
CN113237492A (en) * | 2021-04-27 | 2021-08-10 | 北京航天飞行控制中心 | Initial phase stability determination method, storage medium and electronic device |
CN117494908A (en) * | 2023-12-29 | 2024-02-02 | 宁波港信息通信有限公司 | Port cargo throughput prediction method and system based on big data |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH05216195A (en) * | 1992-02-07 | 1993-08-27 | Konica Corp | Pigment image receiving material |
JPH06149863A (en) * | 1992-11-05 | 1994-05-31 | Yukio Tanaka | Time series data analyzer |
-
2008
- 2008-10-31 JP JP2008282192A patent/JP4711356B2/en not_active Expired - Fee Related
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
JPH05216195A (en) * | 1992-02-07 | 1993-08-27 | Konica Corp | Pigment image receiving material |
JPH06149863A (en) * | 1992-11-05 | 1994-05-31 | Yukio Tanaka | Time series data analyzer |
Non-Patent Citations (2)
Title |
---|
JPN6010073797, 三宅浩次監修, 生物リズムの構造 —MemCalcによる生物時系列データの解析—, 19920115, 初版, pp.19−39, 富士書院 * |
JPN6010073800, 常磐野和男、外2名, 最大エントロピー法による時系列解析 —MemCalcの理論と実際—, 20020625, 第1刷, pp.79−99, 北海道大学図書刊行会 * |
Cited By (7)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN107219840A (en) * | 2017-05-05 | 2017-09-29 | 浙江理工大学 | Towards the regulating valve nonlinear characteristic detection method and system of Natural Gas Station |
JP2019191630A (en) * | 2018-04-18 | 2019-10-31 | 国立大学法人京都大学 | Analysis device, analysis method and computer program |
JP7088538B2 (en) | 2018-04-18 | 2022-06-21 | 国立大学法人京都大学 | Analyst, analysis method, and computer program |
CN113237492A (en) * | 2021-04-27 | 2021-08-10 | 北京航天飞行控制中心 | Initial phase stability determination method, storage medium and electronic device |
CN113237492B (en) * | 2021-04-27 | 2022-05-31 | 北京航天飞行控制中心 | Initial phase stability determination method, storage medium and electronic device |
CN117494908A (en) * | 2023-12-29 | 2024-02-02 | 宁波港信息通信有限公司 | Port cargo throughput prediction method and system based on big data |
CN117494908B (en) * | 2023-12-29 | 2024-03-22 | 宁波港信息通信有限公司 | Port cargo throughput prediction method and system based on big data |
Also Published As
Publication number | Publication date |
---|---|
JP4711356B2 (en) | 2011-06-29 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
Li et al. | Short-term load-forecasting method based on wavelet decomposition with second-order gray neural network model combined with ADF test | |
KR101848958B1 (en) | Improving tool performance by linking spectroscopic information with tool operational parameters and material measurements | |
Satapathy et al. | Empirical assessment of machine learning models for agile software development effort estimation using story points | |
JP4711356B2 (en) | Two-series time-series data analysis apparatus and computer-readable recording medium having recorded two-series time-series data analysis program | |
Wang et al. | Smoothing splines with varying smoothing parameter | |
US20110238606A1 (en) | Kernel regression system, method, and program | |
WO2008151098A1 (en) | Simulating machine and method for determining sensitivity of a system output to changes in underlying system parameters | |
CN112884193A (en) | Prediction device, prediction method, and recording medium | |
Sharma et al. | An approximate analysis of quasi-periodic systems via Floquét theory | |
JPWO2020111258A1 (en) | Virtual measuring device, virtual measuring method and virtual measuring program | |
Cheng et al. | Source contribution evaluation of mechanical vibration signals via enhanced independent component analysis | |
US20120323985A1 (en) | Alignment of multiple liquid chromatography-mass spectrometry runs | |
KR20210092575A (en) | Semiconductor device for compressing a neural network based on a target performance | |
JP4711355B2 (en) | Time series data analyzer | |
CN104345049B (en) | Threshold correction method applied to wavelet threshold noise reduction of laser-induced breakdown spectroscopy | |
Reutskiy | The method of external excitation for solving generalized Sturm–Liouville problems | |
JP2017198620A (en) | Abnormality diagnosis device and abnormality diagnosis method | |
JP2015004573A (en) | Frequency analyzer | |
Shardt | Data quality assessment for closed-loop system identification and forecasting with application to soft sensors | |
EP3768895B1 (en) | Asphalt density estimation system, and related method of reducing signal noise | |
Zhang et al. | Optimal compression of a generalized Prandtl-Ishlinskii operator in hysteresis modeling | |
JP2010185682A (en) | General harmonic analyzer and frequency analyzer | |
Thuan et al. | Efficient bearing fault diagnosis with neural network search and parameter quantization based on vibration and temperature | |
US6600566B1 (en) | High-order high-frequency rough surface scattering solver | |
RU2595929C2 (en) | Method and apparatus for compressing data depending on time signal |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
A621 | Written request for application examination |
Free format text: JAPANESE INTERMEDIATE CODE: A621 Effective date: 20101022 |
|
A871 | Explanation of circumstances concerning accelerated examination |
Free format text: JAPANESE INTERMEDIATE CODE: A871 Effective date: 20101105 |
|
A975 | Report on accelerated examination |
Free format text: JAPANESE INTERMEDIATE CODE: A971005 Effective date: 20101208 |
|
A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20101227 |
|
A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20110225 |
|
A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A821 Effective date: 20110225 |
|
A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20110317 |
|
A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20110317 |
|
R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
LAPS | Cancellation because of no payment of annual fees |