WO2014141344A1 - データ予測装置 - Google Patents

データ予測装置 Download PDF

Info

Publication number
WO2014141344A1
WO2014141344A1 PCT/JP2013/007424 JP2013007424W WO2014141344A1 WO 2014141344 A1 WO2014141344 A1 WO 2014141344A1 JP 2013007424 W JP2013007424 W JP 2013007424W WO 2014141344 A1 WO2014141344 A1 WO 2014141344A1
Authority
WO
WIPO (PCT)
Prior art keywords
state model
model
steady
series data
time
Prior art date
Application number
PCT/JP2013/007424
Other languages
English (en)
French (fr)
Inventor
裕志 吉田
Original Assignee
日本電気株式会社
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by 日本電気株式会社 filed Critical 日本電気株式会社
Priority to JP2015505088A priority Critical patent/JP6337881B2/ja
Priority to US14/775,485 priority patent/US20160042101A1/en
Publication of WO2014141344A1 publication Critical patent/WO2014141344A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation

Definitions

  • the present invention relates to a data prediction apparatus, and more particularly to a data prediction apparatus that predicts time-series data values.
  • communication networks such as the Internet network and mobile packet network
  • communication services are provided mainly on a best-effort basis, but communication throughput, which is the data size (data amount) delivered (transmitted) per unit time, is severely affected by cross traffic and radio wave conditions. Can vary. For this reason, for example, on the service provider side, it is necessary to take a countermeasure in advance by predicting the communication throughput, and a communication throughput prediction apparatus for predicting such communication throughput has been developed.
  • model parameters of a mathematical model are determined from past time series data, and a predicted value is calculated based on the mathematical model.
  • Non-Patent Document 1 As another communication throughput prediction apparatus, there is a communication throughput prediction apparatus described in Non-Patent Document 1.
  • a fluctuation process steady process / unsteady process
  • a mixed model in which a steady process model and an unsteady process model are mixed is constructed based on the discrimination history.
  • a probability distribution probability density function
  • a probability spread probability spread
  • the communication throughput related to communication according to TCP / IP has various factors (for example, end-to-end delay, packet loss, cross traffic, radio wave intensity in wireless communication, etc.). It varies from moment to moment due to the complex action of.
  • model parameters of a mathematical model are determined from past time series data, and a predicted value is calculated based on the mathematical model.
  • the communication throughput fluctuation process stationary process / non-stationary process
  • the probability distribution probability density function
  • any of the above prediction techniques uses a time series model described by a recurrence formula (difference equation) as a prediction model. For this reason, there is a problem that a prediction model cannot be accurately constructed unless the time interval of each point of time series data of past communication throughput observed is equal. Therefore, when the past time-series data of communication throughput is at unequal intervals, the future communication throughput cannot be accurately predicted. Such a problem occurs not only in the prediction of communication throughput, but also in the case of predicting values of all time series data.
  • an object of the present invention is to solve the above-described problem that the value of time series data cannot be predicted with high accuracy.
  • a data prediction apparatus Data observation means for observing time-series data values; A steady-state model representing the time-series data when the fluctuation process of the time-series data is a stationary process, and a non-steady-state model representing the time-series data when the fluctuation process of the time-series data is an unsteady process.
  • Model identification means for identifying each with a stochastic differential equation model based on observed past time series data, Likelihood calculation means for calculating likelihoods that are values representing the likelihood of the steady-state model and the non-steady-state model based on observed past time-series data, A mixing ratio calculating means for calculating a mixing ratio of the steady state model and the unsteady state model based on the likelihood of each of the steady state model and the unsteady state model; A probability distribution predicting means for predicting a probability distribution of time series data based on a prediction model obtained by mixing the steady state model and the non-steady state model according to the mixture ratio;
  • the configuration is as follows.
  • the program which is the other form of this invention is: In the information processing device, Data observation means for observing time-series data values; A steady-state model representing the time-series data when the fluctuation process of the time-series data is a stationary process, and a non-steady-state model representing the time-series data when the fluctuation process of the time-series data is an unsteady process.
  • Model identification means for identifying each with a stochastic differential equation model based on observed past time series data
  • Likelihood calculation means for calculating likelihoods that are values representing the likelihood of the steady-state model and the non-steady-state model based on observed past time-series data
  • a mixing ratio calculating means for calculating a mixing ratio of the steady state model and the unsteady state model based on the likelihood of each of the steady state model and the unsteady state model
  • a probability distribution predicting means for predicting a probability distribution of time series data based on a prediction model obtained by mixing the steady state model and the non-steady state model according to the mixture ratio; It is a program for realizing.
  • a data prediction method includes: Observe the value of time series data, A steady-state model representing the time-series data when the fluctuation process of the time-series data is a stationary process, and a non-steady-state model representing the time-series data when the fluctuation process of the time-series data is an unsteady process. , Identify each with a stochastic differential equation model based on the observed historical time series data, A likelihood that is a value representing the likelihood of the steady-state model and the non-steady-state model is calculated based on observed past time-series data, respectively.
  • a mixture ratio between the steady-state model and the non-steady-state model is calculated, Predicting a probability distribution of time-series data based on a prediction model obtained by mixing the steady-state model and the non-steady-state model according to the mixing ratio;
  • the present invention is configured as described above, so that the value of time series data can be predicted with high accuracy.
  • FIG. 1 is a functional block diagram showing the configuration of the data prediction apparatus.
  • FIG. 2 is a graph showing information used in the data prediction apparatus.
  • FIG. 3 is a schematic diagram showing a probability distribution of data to be predicted.
  • FIG. 4 is a graph comparing the data prediction accuracy in the present embodiment with the data prediction accuracy in other techniques.
  • the data prediction device 1 in the present invention is a general information processing device including an arithmetic device and a storage device.
  • the data predicting apparatus 1 is constructed by incorporating a program into a computing device, and is constructed by a data observing unit 11, a steady stochastic differential equation model identifying unit 12, a non-stationary stochastic differential equation model identifying unit. 13, a likelihood calculation unit 14, a likelihood ratio test unit 15, a mixture ratio calculation unit 16, and a probability distribution prediction unit 17.
  • a data observing unit 11 a steady stochastic differential equation model identifying unit 12
  • a non-stationary stochastic differential equation model identifying unit. 13 a likelihood calculation unit 14
  • a likelihood ratio test unit a mixture ratio calculation unit 16
  • a probability distribution prediction unit 17 a probability distribution prediction unit
  • the data observation unit 11 observes the target time series data ⁇ x t ⁇ .
  • the time series data is an observed data string of random variables that change with time.
  • the time series data that is the target in the data prediction apparatus is not limited to the communication throughput, and may be any time series data.
  • the time interval between adjacent data of the observed time series data needs to be equal.
  • the data prediction apparatus according to the present invention may have unequal time intervals between adjacent data as in the above example. This is because, as will be described later, a model of data at a predetermined time is identified by a stochastic differential equation model.
  • the stationary stochastic differential equation model identification unit 12 (model identification unit), based on the time series data observed by the data observation unit 11 described above, obtains time series data when the variation process of the time series data is a stationary process. Identify the stochastic differential equation model (stationary stochastic differential equation model (steady state model)) to represent.
  • Equation 1 a stochastic differential equation model described by Equation 1 is used as a stochastic differential equation model representing time series data.
  • Equation 1 is a stochastic differential equation model in which the difference is replaced with a differential with respect to the time series model described by the recurrence formula (difference equation) of Non-Patent Document 1 described above. In this way, by making the time interval of the time series model close to infinity, a more accurate data prediction value can be obtained even if the observed time series data is unequal.
  • the stochastic differential equation model expressed by Equation 1 is a stationary process when the real constant a is a> 0, and is a non-stationary process when a ⁇ 0.
  • the stationary stochastic differential equation model identifying unit 12 identifies a stationary stochastic differential equation model of a> 0 in Equation 1. This is equivalent to estimating a, b, and ⁇ which are parameters of the steady stochastic differential equation model of Equation 1.
  • a method for identifying a stationary stochastic differential equation model will be described in detail.
  • Equation 1 The stochastic differential equation model expressed by Equation 1 is a stochastic process called the Orstein-Uhlenbeck process.
  • a, b, and ⁇ are constants, it is called a Vasikek model, and a general solution is obtained.
  • x s is observed at time s
  • Equation 2 the general solution of x t at the subsequent time t (> s) is expressed by Equation 2.
  • Equation 3 the conditional expected value and conditional variance of x t at the subsequent time t (> s) are calculated by Equation 3 and Equation 4, respectively.
  • Equation 2 Since the Orstein-Uhlenbeck process is a class of the Gaussian process, the probability distribution at each time of the general solution expressed by Equation 2 is a Gaussian distribution. Therefore, if E [x t
  • the stationary stochastic differential equation model identification unit 12 aims to estimate model parameters a, b, and ⁇ .
  • a method for estimating the model parameters a, b, and ⁇ using a maximum likelihood estimation method will be described.
  • the likelihood function L is also a, b, It is a function of ⁇ .
  • a, b, and ⁇ that maximize the likelihood function L are obtained.
  • Equation 7 is obtained.
  • ⁇ t i t i ⁇ t i ⁇ 1 .
  • Maximizing the likelihood function L is equal to maximizing the lnL that is the logarithm of the likelihood function L. Since the first term on the right side of Equation 7 is a term irrelevant to a, b, and ⁇ , the sum of the second term and the third term may be maximized.
  • Equation 8 Equation 9
  • the quasi-Newton method is used as a method for calculating a, b, and ⁇ that minimizes F + G.
  • the specific quasi-Newton processing steps are as follows.
  • ⁇ (F + G) is defined by Equation 11.
  • Step 2 The search step width is obtained in accordance with the Armijo condition shown in Step 2.1 to Step 2.4 below.
  • Step 2.2 If the Armijo condition expressed by Equation 12 is satisfied, go to Step 2.4. Otherwise go to 2.3.
  • Step 3) ⁇ is updated by Equation 13.
  • Step 4 End if the stop condition is satisfied. Otherwise go to step 5.
  • the stop condition there are Expression 14 and Expression 15.
  • Equations 16 and 17 are calculated.
  • Step 6 The matrix B k is updated using Equation 18 (BFGS formula).
  • the Armijo condition is used to determine the step width of the search in step 2, but the Wolfe condition may be used.
  • the BFGS formula matrix B k instead of the BFGS formula matrix B k , an H formula that is calculated based on the inverse matrix H k of Bk may be used.
  • the non-stationary stochastic differential equation model identification unit 13 (model identification means) is a time series when the variation process of the time series data is a non-stationary process based on the time series data observed by the data observation unit 11 described above.
  • a non-stationary stochastic differential equation model (unsteady state model)) that is a stochastic differential equation model representing data is identified. That is, the model parameters of the non-stationary stochastic differential equation model are estimated.
  • the stochastic differential equation that is the base of the time-series data model is Equation 1, and this stochastic differential equation is non-stationary when a ⁇ 0.
  • this stochastic differential equation is non-stationary when a ⁇ 0.
  • the stochastic differential equation model of Equation 19 is equivalent to the Brownian motion model, and there is only one model parameter ⁇ . Therefore, in order to identify the unsteady stochastic differential equation model, ⁇ may be estimated.
  • ⁇ is estimated using the maximum likelihood estimation method.
  • the general solution of the unsteady stochastic differential equation model of Equation 19 is Equation 20.
  • conditional expectation, conditional variance, and conditional probability density function of x t at time t (> s) after x s is observed at time s are as shown in Equations 21, 22, and 23.
  • Equation 25 ⁇ that maximizes the logarithm lnL of the likelihood function L of Equation 24 is calculated.
  • the ⁇ is obtained analytically and is given by Equation 25.
  • the likelihood calculating unit 14 is a likelihood that is a value representing the likelihood of each stochastic differential equation model identified by the stationary stochastic differential equation model identifying unit 12 and the non-stationary stochastic differential equation model identifying unit 13. The degree is calculated based on the observed time series data.
  • the likelihood of the stationary stochastic differential equation model can be obtained by calculating based on Equation 6 and the likelihood of the non-stationary stochastic differential equation model based on Equation 24, respectively.
  • the likelihood ratio test unit 15 (test unit) is observed based on the ratio between the likelihood of the stationary stochastic differential equation model calculated by the likelihood calculating unit 14 and the likelihood of the non-stationary stochastic differential equation model.
  • the hypothesis test is performed to determine whether the time series data fits the stationary stochastic differential equation model or the non-stationary stochastic differential equation model.
  • the hypothesis that “the observed time series data is data generated from a non-stationary stochastic differential equation model” is tested in the null hypothesis.
  • the alternative hypothesis is that the observed time series data is data generated from a stationary stochastic differential equation model.
  • R (Equation 27) obtained by multiplying the logarithm of the likelihood ratio ⁇ (Equation 26) defined below by ⁇ 2 is used for the test, where L s is a stationary stochastic differential equation model. (Sup.6), and sup ⁇ L s ⁇ is the upper limit.
  • L n is the likelihood (Expression 24) of the non-stationary stochastic differential equation model, and sup ⁇ L n ⁇ is the upper limit.
  • the likelihood calculated by the likelihood ratio test unit 15 may be used for sup ⁇ L s ⁇ and sup ⁇ L n ⁇ , respectively. This is because the likelihood calculated by the likelihood ratio test unit 15 is a likelihood calculated based on a model parameter that maximizes each likelihood function (Equation 6 and Equation 24), and the likelihood is considered to be an upper limit. Because it is good.
  • the likelihood upper limit sup ⁇ L s ⁇ of the stationary stochastic differential equation model is always greater than or equal to the upper limit sup ⁇ L n ⁇ of the likelihood of the non-stationary stochastic differential equation model (sup ⁇ L s ⁇ ⁇ sup ⁇ L n). ⁇ ). This is because the stationary stochastic differential equation model has three model parameters (a, b, and ⁇ ), whereas the non-stationary stochastic differential equation model has one model parameter ( ⁇ only). Therefore, as in Equation 28, the statistic R is a non-negative real number.
  • the likelihood sup ⁇ L s ⁇ of the stationary stochastic differential equation model is the likelihood of the nonstationary stochastic differential equation model. becomes larger than the degree sup ⁇ L n ⁇ , as a result, statistics values of R by utilizing the fact that larger, statistic R if is becomes greater than a predetermined value, the alternative hypothesis was rejected the null hypothesis (Hypothesis of stationary stochastic differential equation model) is adopted. On the other hand, if the statistic R falls below a predetermined value, the null hypothesis is accepted without being rejected.
  • Threshold of whether to reject the null hypothesis is determined by the distribution of the statistic R when the null hypothesis is correct (this is called the null distribution) and a predetermined significance level. Since it is difficult to obtain the null distribution analytically, in this embodiment, the distribution obtained by the Monte Carlo simulation is used.
  • FIG. 2 shows a null distribution (cumulative distribution function) obtained by Monte Carlo simulation.
  • the null distribution is a distribution obtained by repeating the trial of generating 100 points of time series data and calculating the statistic R under the null hypothesis (non-stationary stochastic differential equation model) 3 million times. .
  • the null hypothesis can be rejected with R> 7.6 when the significance level is 0.1, R> 9.2 when the significance level is 0.05, and R> 12.8 when the significance level is 0.01.
  • the likelihood ratio test unit 15 prepares in advance a threshold value obtained based on the null distribution and significance level, or the null distribution and significance level (for example, when the significance level is 0.1, the threshold is 7.6). Then, the statistic R is calculated from the observed time-series data based on the formulas 26 and 27, and the hypothesis that the model is a stationary stochastic differential equation model is adopted based on the statistic R and the threshold value. Accept the assumption that it is a stochastic differential equation model.
  • equation 30 adopts an exponential load moving average lambda t of the u t to the mixing ratio.
  • ⁇ is a smoothing coefficient for exponential load movement, and 0 ⁇ ⁇ ⁇ 1.
  • the stationary stochastic differential equation model and the non-stationary stochastic differential equation model are mixed. From the definition of Equation 29, the ratio of the non-stationary stochastic differential equation model matches ⁇ t .
  • the probability distribution predicting unit 17 includes the mixing ratio calculated by the mixing ratio calculating unit 16, the steady stochastic differential equation model identified by the steady stochastic differential equation model identifying unit 12 based on the mixing ratio, The probability distribution of future data is predicted from the unsteady stochastic differential equation model identified by the unsteady stochastic differential equation model identifying unit 13.
  • Equation 31 The probability density function of the random variable in the steady stochastic differential equation model expressed by Equation 5 is replaced with f (x t ), and the probability density function of the random variable in the non-stationary stochastic differential equation model expressed by Equation 23 is changed to g.
  • (x t ) the probability density function of the random variable in the non-stationary stochastic differential equation model expressed by Equation 23 is changed to g.
  • Equation 31 is a mixed normal distribution in which two normal distributions are mixed, and the expected value E mix [x t ] and the variance V mix [x t ] are calculated as in Equations 32 and 33.
  • E s [x t ] and V s [x t ] are the expected value and variance of x t in the stationary stochastic differential equation model
  • E n [x t ] and V n [x t ] are non-stationary.
  • the stochastic diffusion represented by Equation 34 is a value obtained by adding or subtracting a constant multiple ( ⁇ times) of the standard deviation from the expected value.
  • FIG. 3 is a schematic diagram showing the probability density function, expected value, and stochastic diffusion of the prediction model. Stochastic diffusion spreads over time, which represents the uncertainty of the predicted value of the data over time. In the stochastic diffusion, the spread increases as the ratio of the non-stationary stochastic differential equation model increases, and the spread decreases as the ratio of the stationary stochastic differential equation model increases.
  • the prediction accuracy is shown in FIG.
  • a diffusion value was obtained from a histogram of variations in actual data values, and a value obtained by subtracting an error (%) from the predicted stochastic diffusion from 100% was taken as a predicted value.
  • the data to be predicted is time-series data of communication throughput in the mobile network, and the time interval of each data is unequal-interval time-series data that follows an exponential distribution of 2 seconds on average. It can be seen that the prediction method using the stochastic differential equation model has higher prediction accuracy.
  • Data observation means 101 for observing time-series data values Data observation means 101 for observing time-series data values; A steady-state model representing the time-series data when the fluctuation process of the time-series data is a stationary process, and a non-steady-state model representing the time-series data when the fluctuation process of the time-series data is an unsteady process.
  • Model identifying means 102 for identifying each with a stochastic differential equation model based on observed past time series data
  • Likelihood calculating means 103 for calculating likelihoods that are values representing the likelihood of the steady-state model and the non-steady-state model based on observed past time-series data
  • a mixing ratio calculating means 104 for calculating a mixing ratio of the steady state model and the unsteady state model based on the likelihood of each of the steady state model and the unsteady state model
  • a probability distribution prediction unit 105 that predicts a probability distribution of time-series data based on a prediction model obtained by mixing the steady-state model and the non-steady-state model according to the mixture ratio
  • a data prediction apparatus 100 comprising:
  • Appendix 2 A data prediction apparatus according to appendix 1, wherein The model identifying means identifies the steady state model and the non-steady state model with different stochastic differential equation models, Data prediction device.
  • Appendix 3 A data prediction apparatus according to appendix 1 or 2, The model identifying means identifies the steady state model with a Vasikek model and identifies the steady state model with a Brownian motion model; Data prediction device.
  • Appendix 4 A data prediction apparatus according to any one of appendices 1 to 3, Based on the ratio between the likelihood of the steady-state model and the likelihood of the non-steady-state model, whether the observed time-series data matches the steady-state model or the non-steady-state model It has a verification means to verify, The mixture ratio calculating means calculates the mixture ratio of the steady state model and the unsteady state model based on the result of the test; Data prediction device.
  • Appendix 5 A data prediction apparatus according to appendix 4, wherein The test means performs a hypothesis test with the null hypothesis that the observed time series data fits the non-steady state model and the alternative hypothesis that the observed time series data fits the steady state model. Do, Data prediction device.
  • the mixing ratio calculation means sets a variable that becomes “0” when the result of the test matches the steady state model and becomes “1” when the non-steady state model matches, and smoothes the variable Calculated as the mixing ratio, Data prediction device.
  • Data observation means for observing time-series data values; A steady-state model representing the time-series data when the fluctuation process of the time-series data is a stationary process, and a non-steady-state model representing the time-series data when the fluctuation process of the time-series data is an unsteady process.
  • Model identification means for identifying each with a stochastic differential equation model based on observed past time series data, Likelihood calculation means for calculating likelihoods that are values representing the likelihood of the steady-state model and the non-steady-state model based on observed past time-series data, A mixing ratio calculating means for calculating a mixing ratio of the steady state model and the unsteady state model based on the likelihood of each of the steady state model and the unsteady state model; A probability distribution predicting means for predicting a probability distribution of time series data based on a prediction model obtained by mixing the steady state model and the non-steady state model according to the mixture ratio; A program to realize
  • the model identifying means identifies the steady state model with a Vasikek model and identifies the steady state model with a Brownian motion model; program.
  • a mixture ratio between the steady-state model and the non-steady-state model is calculated, Predicting a probability distribution of time-series data based on a prediction model obtained by mixing the steady-state model and the non-steady-state model according to the mixing ratio; Data prediction method device.
  • the above-described program is stored in a storage device or recorded on a computer-readable recording medium.
  • the recording medium is a portable medium such as a flexible disk, an optical disk, a magneto-optical disk, and a semiconductor memory.

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Pure & Applied Mathematics (AREA)
  • Mathematical Optimization (AREA)
  • General Engineering & Computer Science (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Databases & Information Systems (AREA)
  • Evolutionary Biology (AREA)
  • Software Systems (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Operations Research (AREA)
  • Probability & Statistics with Applications (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Algebra (AREA)
  • Computer Hardware Design (AREA)
  • Geometry (AREA)
  • Evolutionary Computation (AREA)
  • Complex Calculations (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

　本発明のデータ予測装置は、時系列データの値を観測するデータ観測手段と、定常状態モデルと非定常状態モデルとを観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定するモデル同定手段と、定常状態モデルと非定常状態モデルとの尤もらしさを表す値である尤度をそれぞれ算出する尤度算出手段と、定常状態モデルと非定常状態モデルとのそれぞれの尤度に基づいて定常状態モデルと非定常状態モデルとの混合比を算出する混合比算出手段と、混合比に従って定常状態モデルと非定常状態モデルとを混合して得られる予測モデルに基づいて時系列データの確率分布を予測する確率分布予測手段と、を備える。

Description

データ予測装置
　本発明は、データ予測装置にかかり、特に、時系列データの値を予測するデータ予測装置に関する。
　クラウドサービスの普及により、インターネット網やモバイルパケット網などの通信ネットワークを介した通信量が増加している。そして、通信ネットワークでは、主にベストエフォート型で通信サービスが提供されるが、クロストラヒックや電波状態によって、単位時間あたりに配信（伝送）されるデータサイズ（データの量）である通信スループットを激しく変動しうる。このため、例えば、サービス事業者側では、通信スループットを予測して事前に対策を行う必要があり、このような通信スループットを予測する通信スループット予測装置が開発されている。
　この種の通信スループット予測装置の１つとして、特許文献１に記載されている予測装置がある。この予測装置では、過去の時系列データから数学モデル（線形・非線形混合モデル）のモデルパラメータを決定し、当該数学モデルに基づいて予測値を算出する。
　また、別の通信スループット予測装置の１つとして、非特許文献１に記載されている通信スループット予測装置がある。この予測装置では、通信スループットの変動過程（定常過程／非定常過程）を判別し、かかる判別履歴に基づいて、定常過程モデルと非定常過程モデルを混合した混合モデルを構築し、当該混合モデルに基づいて未来の通信スループットの確率分布（確率密度関数）を算出し、かかる確率密度関数から未来の通信スループットの確率的な広がり（確率的拡散）を算出する。
特開２０１２－１２２８５号公報
吉田裕志，里田浩三，「アプリケーションレベルでのＴＣＰスループットの定常性解析と予測モデル構築」，信学技報，ｖｏｌ．１１２，ｎｏ．３５２，ＩＮ２０１２－１２８，ｐｐ．３９－４４，２０１２年１２月．
　ところで、ＴＣＰ／ＩＰ（Ｔｒａｎｓｍｉｓｓｉｏｎ　Ｃｏｎｔｒｏｌ　Ｐｒｏｔｏｃｏｌ／Ｉｎｔｅｒｎｅｔ　Ｐｒｏｔｏｃｏｌ）に従った通信に係る通信スループットは、様々な要因（例えば，Ｅｎｄ－ｔｏ－Ｅｎｄ遅延、パケットロス、クロストラヒック、及び無線通信における電波強度等）が複雑に作用することによって、時々刻々変動する。
　このような状況に対して、上記特許文献１では、過去の時系列データから数学モデル（線形・非線形混合モデル）のモデルパラメータを決定し、当該数学モデルに基づいて予測値を算出する。また、上記非特許文献１では、上述のように時々刻々変動する通信スループットの変動過程（定常過程／非定常過程）を、観測された過去の通信スループットの時系列データに基づいて判別し、観測された過去の通信スループットの時系列データ及び判別履歴に基づいて定常過程モデルと非定常過程モデルを混合した混合モデルを構築し、当該混合モデルに基づいて未来の通信スループットの確率分布（確率密度関数）を予測することができる。
　しかしながら、上記のいずれの予測技術も、予測モデルとして漸化式（差分方程式）で記述された時系列モデルを用いている。このため、観測された過去の通信スループットの時系列データの各点の時間間隔が等間隔でなければ、正確に予測モデルを構築できない、という問題がある。従って、通信スループットの過去の時系列データが不等間隔であった場合には、正確に未来の通信スループットを予測することができない。また、かかる問題は、通信スループットの予測に限らず、あらゆる時系列データの値を予測する場合も同様に生じる。
　このため、本発明の目的は、上述した課題である、時系列データの値の予測を高精度に行うことができない、という問題を解決することにある。
　本発明の一形態であるデータ予測装置は、
　時系列データの値を観測するデータ観測手段と、
　時系列データの変動過程が定常過程である場合における当該時系列データを表す定常状態モデルと、時系列データの変動過程が非定常過程である場合における当該時系列データを表す非定常状態モデルとを、観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定するモデル同定手段と、
　前記定常状態モデルと前記非定常状態モデルとの尤もらしさを表す値である尤度を、観測された過去の時系列データに基づいてそれぞれ算出する尤度算出手段と、
　前記定常状態モデルと前記非定常状態モデルとのそれぞれの前記尤度に基づいて、前記定常状態モデルと前記非定常状態モデルとの混合比を算出する混合比算出手段と、
　前記混合比に従って前記定常状態モデルと前記非定常状態モデルとを混合して得られる予測モデルに基づいて、時系列データの確率分布を予測する確率分布予測手段と、
を備えた、
という構成をとる。
　また、本発明の他の形態であるプログラムは、
　情報処理装置に、
　時系列データの値を観測するデータ観測手段と、
　時系列データの変動過程が定常過程である場合における当該時系列データを表す定常状態モデルと、時系列データの変動過程が非定常過程である場合における当該時系列データを表す非定常状態モデルとを、観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定するモデル同定手段と、
　前記定常状態モデルと前記非定常状態モデルとの尤もらしさを表す値である尤度を、観測された過去の時系列データに基づいてそれぞれ算出する尤度算出手段と、
　前記定常状態モデルと前記非定常状態モデルとのそれぞれの前記尤度に基づいて、前記定常状態モデルと前記非定常状態モデルとの混合比を算出する混合比算出手段と、
　前記混合比に従って前記定常状態モデルと前記非定常状態モデルとを混合して得られる予測モデルに基づいて、時系列データの確率分布を予測する確率分布予測手段と、
を実現させるためのプログラムである。
　また、本発明の他の形態であるデータ予測方法は、
　時系列データの値を観測し、
　時系列データの変動過程が定常過程である場合における当該時系列データを表す定常状態モデルと、時系列データの変動過程が非定常過程である場合における当該時系列データを表す非定常状態モデルとを、観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定し、
　前記定常状態モデルと前記非定常状態モデルとの尤もらしさを表す値である尤度を、観測された過去の時系列データに基づいてそれぞれ算出し、
　前記定常状態モデルと前記非定常状態モデルとのそれぞれの前記尤度に基づいて、前記定常状態モデルと前記非定常状態モデルとの混合比を算出し、
　前記混合比に従って前記定常状態モデルと前記非定常状態モデルとを混合して得られる予測モデルに基づいて、時系列データの確率分布を予測する、
という構成をとる。
　本発明は、以上のように構成されることにより、時系列データの値を高精度に予測することができる。
本発明の実施形態１におけるデータ予測装置の構成を示す機能ブロック図である。 図１に開示した尤度比検定部において仮説検定に用いる帰無分布（累積分布関数）のグラフである。 図１に開示したデータ予測装置で予測した未来のデータの確率分布の模式図である。 本発明の実施形態１におけるデータ予測装置におけるデータ予測精度と、他の技術におけるデータ予測精度と、を比較したグラフである。 本発明の付記１におけるデータ予測装置の構成を示すブロック図である。
　＜実施形態１＞
　本発明の第１の実施形態を、図１乃至図４を参照して説明する。図１は、データ予測装置の構成を示す機能ブロック図である。図２は、データ予測装置で使用する情報を示すグラフである。図３は、予測するデータの確率分布を示す模式図である。図４は、本実施形態におけるデータ予測精度と他の技術におけるデータ予測精度とを比較したグラフである。
　本発明におけるデータ予測装置１は、演算装置と記憶装置とを備えた一般的な情報処理装置である。そして、データ予測装置１は、図１に示すように、演算装置にプログラムが組み込まれることで構築された、データ観測部１１、定常確率微分方程式モデル同定部１２、非定常確率微分方程式モデル同定部１３、尤度算出部１４、尤度比検定部１５、混合比算出部１６、確率分布予測部１７、を備えている。以下、各構成及びその動作について説明する。
　［データ観測部１１］
　データ観測部１１（データ観測手段）は、対象となる時系列データ｛ｘｔ｝を観測する。時系列データとは、時間経過によって変動する確率変数の観測されたデータ列のことである。例えば、対象となる時系列データが通信スループットであって、時刻ｔ＝０［秒］，ｔ＝１．５［秒］、ｔ＝４．１［秒］に、それぞれｘ＝５［Ｍｂｐｓ］、ｘ＝３［Ｍｂｐｓ］、ｘ＝７［Ｍｂｐｓ］、という数値が観測された場合、観測された時系列データは、｛ｘ０＝５，ｘ１．５＝３，ｘ４．１＝７｝となる。なお、データ予測装置において対象となる時系列データは、通信スループットであることに限定されず、いかなる時系列データであってもよい。
　ここで、従来のデータ予測装置においては、観測された時系列データの隣接するデータ間の時間間隔が等間隔である必要があった。しかしながら、本発明におけるデータ予測装置は、上述の例のように隣接するデータ間の時間間隔が不等であってもよい。このことは、後述するように、所定の時間におけるデータのモデルを、確率微分方程式モデルで同定していることによる。
　［定常確率微分方程式モデル同定部１２］
　定常確率微分方程式モデル同定部１２（モデル同定手段）は、上述したデータ観測部１１で観測された時系列データに基づいて、当該時系列データの変動過程が定常過程である場合における時系列データを表す確率微分方程式モデル（定常確率微分方程式モデル（定常状態モデル））を同定する。
　ここで、本実施形態では、時系列データを表す確率微分方程式モデルとして、数１で記述される確率微分方程式モデルを用いる。
Figure JPOXMLDOC01-appb-M000001
　上記ｘｔは、対象とする確率変数を示す。また、上記ａ，ｂは実定数、σは正の定数、Ｂｔは標準ブラウン運動、である。数１は、上述した非特許文献１の漸化式（差分方程式）で記述された時系列モデルに対して、差分を微分に置き換えた確率微分方程式モデルである。このように、時系列モデルの時間間隔を無限小に近づけることで、観測された時系列データが不等であっても、より正確なデータ予測値を得ることができる。
　数１で表された確率微分方程式モデルは、実定数ａが、ａ＞０のとき定常過程となり、ａ≦０のとき、非定常過程となることが知られている。このため、定常確率微分方程式モデル同定部１２は、数１においてａ＞０の定常確率微分方程式モデルを同定する。これは、数１の定常確率微分方程式モデルのパラメータであるａ，ｂ，σを推定することに等しい。以下，定常確率微分方程式モデルの同定方法について詳細に説明する。
　数１で表される確率微分方程式モデルは、Ｏｒｎｓｔｅｉｎ－Ｕｈｌｅｎｂｅｃｋ過程と呼ばれる確率過程であり、特に、ａ，ｂ，σが定数であるとき、Ｖａｓｉｃｅｋモデルと呼ばれ、一般解が得られている。時刻ｓでｘｓが観測された時、その後の時刻ｔ（＞ｓ）におけるｘｔの一般解は、数２で表される。
Figure JPOXMLDOC01-appb-M000002
　数２の一般解から、同じく時刻ｓでｘｓが観測された時、その後の時刻ｔ（＞ｓ）におけるｘｔの条件付き期待値及び条件付き分散は、それぞれ数３、数４で計算される。
Figure JPOXMLDOC01-appb-M000003
Figure JPOXMLDOC01-appb-M000004
　Ｏｒｎｓｔｅｉｎ－Ｕｈｌｅｎｂｅｃｋ過程は、ガウス過程のクラスであるから、数２で表される一般解の各時刻における確率分布はガウス分布になる。したがって、数３のＥ［ｘｔ｜ｘｓ］及び数４のＶ［ｘｔ｜ｘｓ］を、改めてμｔ，ｓ及びσ２ ｔ，ｓとおくと、時刻ｓでｘｓが観測された時、その後の時刻ｔ（＞ｓ）におけるｘｔの条件付き確率密度関数は数５で表される。
Figure JPOXMLDOC01-appb-M000005
　上述したとおり、定常確率微分方程式モデル同定部１２は、モデルパラメータであるａ，ｂ，σを推定することを目的とする。本実施形態では、最尤推定法を用いて上記モデルパラメータａ，ｂ，σを推定する方法について説明する。
　まず、ｎ個の過去の時系列データ｛ｘｔ１，ｘｔ２，…，ｘｔｎ｝（ｔ１＜ｔ２＜…＜ｔｎ）が観測されたとする。隣接するデータ間の時間間隔ｔｉ＋１－ｔｉ（ｉ＝１，２，…，ｎ－１）は不等でもよい。定常確率微分方程式モデルの一般解の条件付き確率密度関数は、数５で表されるため、上記ｎ個の過去の時系列データが観測されたときの尤度関数Ｌは、数６のようになる。
Figure JPOXMLDOC01-appb-M000006
　上記数６におけるμｔｉ，ｔｉ－１及びσｔｉ，ｔｉ－１は、それぞれ数３及び数４で表されるとおりａ，ｂ，σの関数であるため、尤度関数Ｌもａ，ｂ，σの関数である。最尤推定法では、この尤度関数Ｌを最大にするａ，ｂ，σを求める。
　しかしながら、解析的には尤度関数Ｌを最大にするａ，ｂ，σを求めることは困難であるため、本実施形態では数値的に尤度関数Ｌを最大にするａ，ｂ，σを求める方法について説明する。
　まず，数６の尤度関数Ｌの対数ｌｎＬを求めると、数７のようになる。ただし、Δｔｉ＝ｔｉ－ｔｉ－１とおいた。
Figure JPOXMLDOC01-appb-M000007
　尤度関数Ｌを最大にすることと、尤度関数Ｌの対数であるｌｎＬを最大にすることとは等しい。数７の右辺第一項は、ａ，ｂ，σに無関係な項であるため、第二項と第三項の和を最大化すればよい。
　ここで，数７の右辺第二項及び第三項それぞれについて、「－１／２」を除いた関数を数８及び数９で定義する。
Figure JPOXMLDOC01-appb-M000008
Figure JPOXMLDOC01-appb-M000009
　結局，尤度関数Ｌを最大することは、上記Ｆ＋Ｇを最小にすることに等しい。Ｆ＋Ｇを最小にするａ，ｂ，σを算出する方法として、本実施形態では準ニュートン法を用いる。具体的な準ニュートンの処理ステップは、以下のとおりである。
（準備）θ＝［ａ　ｂ　σ］Ｔ　（Ｔは転置を表す）とおく。
（ステップ０）適当な初期値θ０を与え，初期Ｂ０は３×３の単位行列とする。
（ステップ１）数１０で表される連立一次方程式を解いて探索方向ベクトルｄを求める。
Figure JPOXMLDOC01-appb-M000010
 ただし、∇（Ｆ＋Ｇ）は数１１で定義される。
Figure JPOXMLDOC01-appb-M000011
（ステップ２）以下のステップ２．１からステップ２．４で示したＡｒｍｉｊｏ条件に従って探索のステップ幅を求める。
（ステップ２．１）βｋ，０＝１，ｉ＝０，０＜ξ＜１，０＜τ＜１とおく。
（ステップ２．２）数１２で表されるＡｒｍｉｊｏ条件を満足するならステップ２．４へ。それ以外は２．３へ。
Figure JPOXMLDOC01-appb-M000012
（ステップ２．３）βｋ，ｉ＋１＝τβｋ，ｉ，ｉ：＝ｉ＋１とおいてステップ２．２へ戻る。
（ステップ２．４）αｋ＝βｋ，ｉとおく。
（ステップ３）数１３でθを更新する。
Figure JPOXMLDOC01-appb-M000013
（ステップ４）停止条件が満足されれば終了。それ以外はステップ５へ。停止条件としては、数１４や数１５がある。
Figure JPOXMLDOC01-appb-M000014
Figure JPOXMLDOC01-appb-M000015
（ステップ５）数１６及び数１７を計算する。
Figure JPOXMLDOC01-appb-M000016
Figure JPOXMLDOC01-appb-M000017
（ステップ６）数１８（ＢＦＧＳ公式）を用いて、行列Ｂｋを更新する。
Figure JPOXMLDOC01-appb-M000018
（ステップ７）ｋ：＝ｋ＋１とおいてステップ１に戻る。
　以上のステップ１からステップ７を実施することで、Ｆ＋Ｇを最小にするθ＝［ａ　ｂ　σ］Ｔを算出することができる．
　上記の準ニュートン法では、ステップ２において探索のステップ幅を求めるのにＡｒｍｉｊｏ条件を用いたが、Ｗｏｌｆｅ条件を用いても良い。また、ＢＦＧＳ公式の行列Ｂｋの代わりに、Ｂｋの逆行列Ｈｋに基づいて計算するＨ公式を用いても良い。
　［非定常確率微分方程式モデル同定部１３］
　非定常確率微分方程式モデル同定部１３（モデル同定手段）は、上述したデータ観測部１１で観測された時系列データに基づいて、当該時系列データの変動過程が非定常過程である場合における時系列データを表す確率微分方程式モデルである非定常確率微分方程式モデル（非定常状態モデル））を同定する。つまり、非定常確率微分方程式モデルのモデルパラメータを推定する。
　ここで、上述したとおり、時系列データのモデルのベースとなる確率微分方程式は、数１であり、この確率微分方程式が非定常となるのは、ａ≦０のときである。しかし、ａ＜０の領域は、急速に無限大に発散する過程となるため、ほとんどの有界な時系列データの予測には不適である。そのため、非定常確率微分方程式モデルとしては、ａ＝０のときのみを考えれば良い。このとき、非定常確率微分方程式モデルは数１９のようになる。
Figure JPOXMLDOC01-appb-M000019
　数１９の確率微分方程式モデルは、ブラウン運動モデルに等しく、モデルパラメータはσの一つだけである。そのため、非定常確率微分方程式モデルを同定するためには、σを推定すればよい。ここでも、定常確率微分方程式モデル同定部１２と同様に、最尤推定法を用いてσを推定する。なお、数１９の非定常確率微分方程式モデルの一般解は、数２０である。
Figure JPOXMLDOC01-appb-M000020
　時刻ｓでｘｓが観測された後の時刻ｔ（＞ｓ）におけるｘｔの条件付き期待、条件付き分散、条件付き確率密度関数は、数２１，２２，２３のようになる。
Figure JPOXMLDOC01-appb-M000021
Figure JPOXMLDOC01-appb-M000022
Figure JPOXMLDOC01-appb-M000023
　このとき，ｎ個の過去の時系列データ｛ｘｔ１，ｘｔ２，…，ｘｔｎ｝（ｔ１＜ｔ２＜…＜ｔｎ）が観測されたときの尤度関数Ｌは、数２４のようになる。ただし、Δｔｉ＝ｔｉ－ｔｉ－１とおいた。
Figure JPOXMLDOC01-appb-M000024
　数２４の尤度関数Ｌの対数ｌｎＬを最大化するσを算出する。当該σは、解析的に求まり、数２５のようになる。
Figure JPOXMLDOC01-appb-M000025
　［尤度算出部１４］
　尤度算出部１４（尤度算出手段）は、上記定常確率微分方程式モデル同定部１２及び上記非定常確率微分方程式モデル同定部１３で同定した各確率微分方程式モデルの尤もらしさを表す値である尤度を、観測された時系列データに基づいてそれぞれ算出する。定常確率微分方程式モデルの尤度は数６、非定常確率微分方程式モデルの尤度は数２４、に基づいてそれぞれ計算することで得られる。
　［尤度比検定部１５］
　尤度比検定部１５（検定手段）は、上記尤度算出部１４で算出した定常確率微分方程式モデルの尤度と、非定常確率微分方程式モデルの尤度と、の比に基づいて、観測された時系列データが、定常確率微分方程式モデルに適合するか非定常確率微分方程式モデルに適合するか仮説検定にかける。
　本実施形態では、「観測された時系列データは、非定常確率微分方程式モデルから生成されたデータである」という仮説を、帰無仮説において検定する。このとき、対立仮説は、「観測された時系列データは定常確率微分方程式モデルから生成されたデータである」となる。
　具体的に、本実施形態では、以下で定義する尤度比Λ（数２６）の対数に－２をかけたＲ（数２７）を検定に用いる、ただし、Ｌｓは、定常確率微分方程式モデルの尤度（数６）であり、ｓｕｐ｛Ｌｓ｝は、その上限である。また、Ｌｎは、非定常確率微分方程式モデルの尤度（数２４）であり、ｓｕｐ｛Ｌｎ｝は、その上限である。
Figure JPOXMLDOC01-appb-M000026
Figure JPOXMLDOC01-appb-M000027
　ｓｕｐ｛Ｌｓ｝及びｓｕｐ｛Ｌｎ｝は、それぞれ尤度比検定部１５で算出した尤度を用いればよい。なぜなら、尤度比検定部１５で算出した尤度は、各尤度関数（数６と数２４）を最大にするモデルパラメータに基づいて算出される尤度であり、当該尤度は上限と考えてよいからである。
　定常確率微分方程式モデルの尤度の上限ｓｕｐ｛Ｌｓ｝は、必ず非定常確率微分方程式モデルの尤度の上限ｓｕｐ｛Ｌｎ｝以上の値となる（ｓｕｐ｛Ｌｓ｝≧ｓｕｐ｛Ｌｎ｝）。これは、定常確率微分方程式モデルのモデルパラメータが３つ（ａとｂとσ）であるのに対し、非定常確率微分方程式モデルのモデルパラメータは１つ（σのみ）だからである。従って、数２８のように、統計量Ｒは非負の実数となる。
Figure JPOXMLDOC01-appb-M000028
　尤度比検定では、帰無仮説（非定常確率微分方程式モデルであるという仮説）が誤りの場合には、定常確率微分方程式モデルの尤度ｓｕｐ｛Ｌｓ｝が非定常確率微分方程式モデルの尤度ｓｕｐ｛Ｌｎ｝に比べて大きくなり、その結果、統計量Ｒの値が大きくなることを利用し、統計量Ｒが所定の値よりも大きくなれば、帰無仮説を棄却して対立仮説（定常確率微分方程式モデルであるという仮説）を採択する。一方、統計量Ｒが所定の値以下になれば、帰無仮説は棄却できずに容認することになる。
　帰無仮説を棄却するか否かの閾値は、帰無仮説が正しい場合の統計量Ｒの分布（これは帰無分布と呼ばれる）と所定の有意水準によって定まる。帰無分布を解析的に求めることは困難であるため、本実施形態では、モンテカルロ・シミュレーションにて求めた分布を用いる。図２に、モンテカルロ・シミュレーションで求めた帰無分布（累積分布関数）を示す。上記帰無分布は、帰無仮説（非定常確率微分方程式モデル）の下で１００点の時系列データを生成して統計量Ｒを算出するという試行を３００万回繰り返して得られた分布である。有意水準０．１の場合はＲ＞７．６、有意水準０．０５の場合はＲ＞９．２、有意水準０．０１の場合はＲ＞１２．８、で帰無仮説を棄却できる。
　尤度比検定部１５は、予め上記帰無分布と有意水準、もしくは帰無分布と有意水準に基づいて得られた閾値（例えば有意水準０．１とき閾値は７．６）を用意しておき、観測された時系列データから数２６，数２７に基づいて統計量Ｒを計算し、当該統計量Ｒと上記閾値に基づいて、定常確率微分方程式モデルであるという仮説を採択するか、非定常確率微分方程式モデルであるという仮定を容認する。
　［混合比算出部１６］
　混合比算出部１６（混合比算出手段）、上記尤度比検定部１５の検定結果の履歴に基づいて、定常確率微分方程式モデル同定部１２で同定した定常確率微分方程式モデルと、上記非定常確率微分方程式モデル同定部１３で同定した非定常確率微分方程式モデルと、を混合する割合を表す混合比を算出する。
　上記尤度比検定部１５にて検定した結果、定常確率微分方程式モデルであることが採択された場合に「０」、非定常確率微分方程式モデルが容認された場合に「１」、をとるような確率変数ｕｔを定義する（数２９）。
Figure JPOXMLDOC01-appb-M000029
　本実施形態では、数３０のように、上記ｕｔの指数荷重移動平均λｔを混合比に採用する。ただし、γは指数荷重移動の平滑化係数であり、０≦γ≦１である。
Figure JPOXMLDOC01-appb-M000030
　得られた混合比λｔに基づいて、定常確率微分方程式モデルと非定常確率微分方程式モデルとを混合する。数２９の定義より、非定常確率微分方程式モデルの割合がλｔに一致する。
　［確率分布予測部１７］
　確率分布予測部１７（確率分布予測手段）は、上記混合比算出部１６で算出した混合比と、当該混合比に基づいて定常確率微分方程式モデル同定部１２で同定した定常確率微分方程式モデルと、非定常確率微分方程式モデル同定部１３で同定した非定常確率微分方程式モデルとから、未来のデータの確率分布を予測する。
　数５で表される定常確率微分方程式モデルにおける確率変数の確率密度関数を改めてｆ（ｘｔ）とおき、数２３で表される非定常確率微分方程式モデルにおける確率変数の確率密度関数を改めてｇ（ｘｔ）とおく。すると、上記混合比算出部１６で算出した混合比λｔに基づいて、混合したモデルにおける確率変数ｘｔの確率密度関数ｈ（ｘｔ）は、数３１で表され、これが未来のデータの確率分布である。
Figure JPOXMLDOC01-appb-M000031
　数３１は二つの正規分布が混合された混合正規分布であり、期待値Ｅｍｉｘ［ｘｔ］及び分散Ｖｍｉｘ［ｘｔ］は、数３２，数３３のように計算される。ただし、Ｅｓ［ｘｔ］及びＶｓ［ｘｔ］は、定常確率微分方程式モデルにおけるｘｔの期待値及び分散であり、Ｅｎ［ｘｔ］及びＶｎ［ｘｔ］は、非定常確率微分方程式モデルにおけるｘｔの期待値及び分散である。
Figure JPOXMLDOC01-appb-M000032
Figure JPOXMLDOC01-appb-M000033
　［発明の効果］
　ここで、未来のデータの値を予測する場合、未来のデータが確率的にどの程度の幅に存在するかという目安が分かれば便利なことがある。この確率的な変動幅を確率的拡散と呼び、数３４で定義する。
Figure JPOXMLDOC01-appb-M000034
　数３４で表される確率的拡散は、期待値から標準偏差の定数倍（α倍）だけ加減した値である。図３は、当該予測モデルの確率密度関数、期待値、確率的拡散を示した模式図である。確率的拡散は、時間経過に従って広がっており、これは時間経過に伴うデータの予測値の不確かさを表す。確率的拡散は、非定常確率微分方程式モデルの割合が高いほど広がりが大きくなり、定常確率微分方程式モデルの割合が高いほど広がりは小さくなる。
　上記確率的拡散の予測精度について、本発明の実施形態の確率微分方程式モデルでの予測方法にて予測した確率的拡散と、従来技術である時系列モデル（漸化式）で予測した確率的拡散との予測精度を、図４に示した。実際のデータ値のばらつきのヒストグラムから拡散値を求め、予測した確率的拡散との誤差（％）を１００％から減じた値を予測値とした。予測対象のデータはモバイルネットワークにおける通信スループットの時系列データであり、各データの時間間隔は平均２秒の指数分布に従う不等間隔時系列データである。確率微分方程式モデルでの予測方法の方が、予測精度が高くなっていることがわかる。
　＜付記＞
　上記実施形態の一部又は全部は、以下の付記のようにも記載されうる。以下、本発明におけるデータ予測装置（図５参照）、プログラム、データ予測方法の構成の概略を説明する。但し、本発明は、以下の構成に限定されない。
（付記１）
　時系列データの値を観測するデータ観測手段１０１と、
　時系列データの変動過程が定常過程である場合における当該時系列データを表す定常状態モデルと、時系列データの変動過程が非定常過程である場合における当該時系列データを表す非定常状態モデルとを、観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定するモデル同定手段１０２と、
　前記定常状態モデルと前記非定常状態モデルとの尤もらしさを表す値である尤度を、観測された過去の時系列データに基づいてそれぞれ算出する尤度算出手段１０３と、
　前記定常状態モデルと前記非定常状態モデルとのそれぞれの前記尤度に基づいて、前記定常状態モデルと前記非定常状態モデルとの混合比を算出する混合比算出手段１０４と、
　前記混合比に従って前記定常状態モデルと前記非定常状態モデルとを混合して得られる予測モデルに基づいて、時系列データの確率分布を予測する確率分布予測手段１０５と、
を備えたデータ予測装置１００。
（付記２）
　付記１に記載のデータ予測装置であって、
　前記モデル同定手段は、前記定常状態モデルと前記非定常状態モデルとをそれぞれ異なる確率微分方程式モデルで同定する、
データ予測装置。
（付記３）
　付記１又は２に記載のデータ予測装置であって、
　前記モデル同定手段は、前記定常状態モデルをＶａｓｉｃｅｋモデルで同定し、前記定常状態モデルをブラウン運動モデルで同定する、
データ予測装置。
（付記４）
　付記１乃至３のいずれかに記載のデータ予測装置であって、
　前記定常状態モデルの前記尤度と前記非定常状態モデルの前記尤度との比に基づいて、観測された時系列データが前記定常状態モデルと前記非定常状態モデルとのいずれに適合するかを検定する検定手段を備え、
　前記混合比算出手段は、前記検定の結果に基づいて前記定常状態モデルと前記非定常状態モデルとの前記混合比を算出する、
データ予測装置。
（付記５）
　付記４に記載のデータ予測装置であって、
　前記検定手段は、観測された時系列データが前記非定常状態モデルに適合することを帰無仮説とし、観測された時系列データが前記定常状態モデルに適合することを対立仮説とする仮説検定を行う、
データ予測装置。
（付記６）
　付記４又は５に記載のデータ予測装置であって、
　前記混合比算出手段は、前記検定の結果、前記定常状態モデルに適合した場合に「０」となり、前記非定常状態モデルに適合した場合に「１」となる変数を設定し、当該変数を平滑化した値を前記混合比として算出する、
データ予測装置。
（付記７）
　情報処理装置に、
　時系列データの値を観測するデータ観測手段と、
　時系列データの変動過程が定常過程である場合における当該時系列データを表す定常状態モデルと、時系列データの変動過程が非定常過程である場合における当該時系列データを表す非定常状態モデルとを、観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定するモデル同定手段と、
　前記定常状態モデルと前記非定常状態モデルとの尤もらしさを表す値である尤度を、観測された過去の時系列データに基づいてそれぞれ算出する尤度算出手段と、
　前記定常状態モデルと前記非定常状態モデルとのそれぞれの前記尤度に基づいて、前記定常状態モデルと前記非定常状態モデルとの混合比を算出する混合比算出手段と、
　前記混合比に従って前記定常状態モデルと前記非定常状態モデルとを混合して得られる予測モデルに基づいて、時系列データの確率分布を予測する確率分布予測手段と、
を実現させるためのプログラム。
（付記８）
　付記７に記載のプログラムであって、
　前記モデル同定手段は、前記定常状態モデルをＶａｓｉｃｅｋモデルで同定し、前記定常状態モデルをブラウン運動モデルで同定する、
プログラム。
（付記９）
　時系列データの値を観測し、
　時系列データの変動過程が定常過程である場合における当該時系列データを表す定常状態モデルと、時系列データの変動過程が非定常過程である場合における当該時系列データを表す非定常状態モデルとを、観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定し、
　前記定常状態モデルと前記非定常状態モデルとの尤もらしさを表す値である尤度を、観測された過去の時系列データに基づいてそれぞれ算出し、
　前記定常状態モデルと前記非定常状態モデルとのそれぞれの前記尤度に基づいて、前記定常状態モデルと前記非定常状態モデルとの混合比を算出し、
　前記混合比に従って前記定常状態モデルと前記非定常状態モデルとを混合して得られる予測モデルに基づいて、時系列データの確率分布を予測する、
データ予測方法装置。
（付記１０）
　付記９に記載のデータ予測方法であって、
　前記定常状態モデルをＶａｓｉｃｅｋモデルで同定し、前記定常状態モデルをブラウン運動モデルで同定する、
データ予測方法。
　なお、上述したプログラムは、記憶装置に記憶されていたり、コンピュータが読み取り可能な記録媒体に記録されている。例えば、記録媒体は、フレキシブルディスク、光ディスク、光磁気ディスク、及び、半導体メモリ等の可搬性を有する媒体である。
　以上、上記実施形態等を参照して本願発明を説明したが、本願発明は、上述した実施形態に限定されるものではない。本願発明の構成や詳細には、本願発明の範囲内で当業者が理解しうる様々な変更をすることができる。
　なお、本発明は、日本国にて２０１３年３月１４日に特許出願された特願２０１３－０５１２０５の特許出願に基づく優先権主張の利益を享受するものであり、当該特許出願に記載された内容は、全て本明細書に含まれるものとする。
１　データ予測装置
１１　データ観測部
１２　定常確率微分方程式モデル同定部
１３　非定常確率微分方程式モデル同定部
１４　尤度算出部
１５　尤度比検定部
１６　混合比算出部
１７　確率分布予測部
１００　データ予測装置
１０１　データ観測手段
１０２　モデル同定手段
１０３　尤度算出手段
１０４　混合比算出手段
１０５　確率分布予測手段
 

Claims (10)

  1. 　時系列データの値を観測するデータ観測手段と、
    　時系列データの変動過程が定常過程である場合における当該時系列データを表す定常状態モデルと、時系列データの変動過程が非定常過程である場合における当該時系列データを表す非定常状態モデルとを、観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定するモデル同定手段と、
    　前記定常状態モデルと前記非定常状態モデルとの尤もらしさを表す値である尤度を、観測された過去の時系列データに基づいてそれぞれ算出する尤度算出手段と、
    　前記定常状態モデルと前記非定常状態モデルとのそれぞれの前記尤度に基づいて、前記定常状態モデルと前記非定常状態モデルとの混合比を算出する混合比算出手段と、
    　前記混合比に従って前記定常状態モデルと前記非定常状態モデルとを混合して得られる予測モデルに基づいて、時系列データの確率分布を予測する確率分布予測手段と、
    を備えたデータ予測装置。
  2. 　請求項１に記載のデータ予測装置であって、
    　前記モデル同定手段は、前記定常状態モデルと前記非定常状態モデルとをそれぞれ異なる確率微分方程式モデルで同定する、
    データ予測装置。
  3. 　請求項１又は２に記載のデータ予測装置であって、
    　前記モデル同定手段は、前記定常状態モデルをＶａｓｉｃｅｋモデルで同定し、前記定常状態モデルをブラウン運動モデルで同定する、
    データ予測装置。
  4. 　請求項１乃至３のいずれかに記載のデータ予測装置であって、
    　前記定常状態モデルの前記尤度と前記非定常状態モデルの前記尤度との比に基づいて、観測された時系列データが前記定常状態モデルと前記非定常状態モデルとのいずれに適合するかを検定する検定手段を備え、
    　前記混合比算出手段は、前記検定の結果に基づいて前記定常状態モデルと前記非定常状態モデルとの前記混合比を算出する、
    データ予測装置。
  5. 　請求項４に記載のデータ予測装置であって、
    　前記検定手段は、観測された時系列データが前記非定常状態モデルに適合することを帰無仮説とし、観測された時系列データが前記定常状態モデルに適合することを対立仮説とする仮説検定を行う、
    データ予測装置。
  6. 　請求項４又は５に記載のデータ予測装置であって、
    　前記混合比算出手段は、前記検定の結果、前記定常状態モデルに適合した場合に「０」となり、前記非定常状態モデルに適合した場合に「１」となる変数を設定し、当該変数を平滑化した値を前記混合比として算出する、
    データ予測装置。
  7. 　情報処理装置に、
    　時系列データの値を観測するデータ観測手段と、
    　時系列データの変動過程が定常過程である場合における当該時系列データを表す定常状態モデルと、時系列データの変動過程が非定常過程である場合における当該時系列データを表す非定常状態モデルとを、観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定するモデル同定手段と、
    　前記定常状態モデルと前記非定常状態モデルとの尤もらしさを表す値である尤度を、観測された過去の時系列データに基づいてそれぞれ算出する尤度算出手段と、
    　前記定常状態モデルと前記非定常状態モデルとのそれぞれの前記尤度に基づいて、前記定常状態モデルと前記非定常状態モデルとの混合比を算出する混合比算出手段と、
    　前記混合比に従って前記定常状態モデルと前記非定常状態モデルとを混合して得られる予測モデルに基づいて、時系列データの確率分布を予測する確率分布予測手段と、
    を実現させるためのプログラム。
  8. 　請求項７に記載のプログラムであって、
    　前記モデル同定手段は、前記定常状態モデルをＶａｓｉｃｅｋモデルで同定し、前記定常状態モデルをブラウン運動モデルで同定する、
    プログラム。
  9. 　時系列データの値を観測し、
    　時系列データの変動過程が定常過程である場合における当該時系列データを表す定常状態モデルと、時系列データの変動過程が非定常過程である場合における当該時系列データを表す非定常状態モデルとを、観測された過去の時系列データに基づいてそれぞれ確率微分方程式モデルで同定し、
    　前記定常状態モデルと前記非定常状態モデルとの尤もらしさを表す値である尤度を、観測された過去の時系列データに基づいてそれぞれ算出し、
    　前記定常状態モデルと前記非定常状態モデルとのそれぞれの前記尤度に基づいて、前記定常状態モデルと前記非定常状態モデルとの混合比を算出し、
    　前記混合比に従って前記定常状態モデルと前記非定常状態モデルとを混合して得られる予測モデルに基づいて、時系列データの確率分布を予測する、
    データ予測方法装置。
  10. 　請求項９に記載のデータ予測方法であって、
    　前記定常状態モデルをＶａｓｉｃｅｋモデルで同定し、前記定常状態モデルをブラウン運動モデルで同定する、
    データ予測方法。
PCT/JP2013/007424 2013-03-14 2013-12-18 データ予測装置 WO2014141344A1 (ja)

Priority Applications (2)

Application Number Priority Date Filing Date Title
JP2015505088A JP6337881B2 (ja) 2013-03-14 2013-12-18 データ予測装置
US14/775,485 US20160042101A1 (en) 2013-03-14 2013-12-18 Data prediction apparatus

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2013-051205 2013-03-14
JP2013051205 2013-03-14

Publications (1)

Publication Number Publication Date
WO2014141344A1 true WO2014141344A1 (ja) 2014-09-18

Family

ID=51536047

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2013/007424 WO2014141344A1 (ja) 2013-03-14 2013-12-18 データ予測装置

Country Status (3)

Country Link
US (1) US20160042101A1 (ja)
JP (1) JP6337881B2 (ja)
WO (1) WO2014141344A1 (ja)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2018139300A1 (ja) * 2017-01-24 2018-08-02 日本電気株式会社 情報処理装置、情報処理方法、及び、情報処理プログラムが記録された記録媒体
WO2019215810A1 (ja) * 2018-05-08 2019-11-14 株式会社日立製作所 データ分析装置、電力潮流解析装置およびデータ分析方法
KR20210092482A (ko) * 2020-01-16 2021-07-26 주식회사 에이젠글로벌 인공지능을 이용한 사기거래 탐지 시스템 및 사기거래 탐지 방법

Families Citing this family (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US9934259B2 (en) 2013-08-15 2018-04-03 Sas Institute Inc. In-memory time series database and processing in a distributed environment
US9892370B2 (en) 2014-06-12 2018-02-13 Sas Institute Inc. Systems and methods for resolving over multiple hierarchies
US9418339B1 (en) * 2015-01-26 2016-08-16 Sas Institute, Inc. Systems and methods for time series analysis techniques utilizing count data sets
CN107665276A (zh) * 2017-09-18 2018-02-06 天津大学 基于符号化模态及转换频次的时间序列复杂性测算方法
US10560313B2 (en) 2018-06-26 2020-02-11 Sas Institute Inc. Pipeline system for time-series data forecasting
US10685283B2 (en) 2018-06-26 2020-06-16 Sas Institute Inc. Demand classification based pipeline system for time-series data forecasting
JP2022072271A (ja) 2020-10-29 2022-05-17 本田技研工業株式会社 情報処理装置、移動体、プログラム、及び情報処理方法
JP7410839B2 (ja) 2020-10-29 2024-01-10 本田技研工業株式会社 情報処理装置、移動体、プログラム、及び情報処理方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
YUJI YOSHIDA ET AL.: "Application of TCP Throughput Prediction to Video Streaming Control and Its Evaluation", IEICE TECHNICAL REPORT, vol. 112, no. 464, 28 February 2013 (2013-02-28), pages 281 - 286 *

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2018139300A1 (ja) * 2017-01-24 2018-08-02 日本電気株式会社 情報処理装置、情報処理方法、及び、情報処理プログラムが記録された記録媒体
JPWO2018139300A1 (ja) * 2017-01-24 2019-11-07 日本電気株式会社 情報処理装置、情報処理方法、及び、情報処理プログラムが記録された記録媒体
WO2019215810A1 (ja) * 2018-05-08 2019-11-14 株式会社日立製作所 データ分析装置、電力潮流解析装置およびデータ分析方法
JPWO2019215810A1 (ja) * 2018-05-08 2021-01-14 株式会社日立製作所 データ分析装置、電力潮流解析装置およびデータ分析方法
KR20210092482A (ko) * 2020-01-16 2021-07-26 주식회사 에이젠글로벌 인공지능을 이용한 사기거래 탐지 시스템 및 사기거래 탐지 방법
KR102412432B1 (ko) * 2020-01-16 2022-06-23 주식회사 에이젠글로벌 인공지능을 이용한 사기거래 탐지 시스템 및 사기거래 탐지 방법

Also Published As

Publication number Publication date
JP6337881B2 (ja) 2018-06-06
US20160042101A1 (en) 2016-02-11
JPWO2014141344A1 (ja) 2017-02-16

Similar Documents

Publication Publication Date Title
JP6337881B2 (ja) データ予測装置
CN110149237A (zh) 一种Hadoop平台计算节点负载预测方法
Ganji et al. Advance first order second moment (AFOSM) method for single reservoir operation reliability analysis: a case study
Azzouz et al. Steady state IBEA assisted by MLP neural networks for expensive multi-objective optimization problems
Dushkin et al. An improved method for predicting the evolution of the characteristic parameters of an information system
Daolio et al. Local optima networks and the performance of iterated local search
CN111460692A (zh) 考虑退化速率相互影响的设备剩余寿命预测方法及系统
US20090138238A1 (en) Sequential fixed-point quantile estimation
Lee et al. Test for parameter change in diffusion processes by cusum statistics based on one-step estimators
Abate et al. Approximate abstractions of stochastic systems: A randomized method
Oreshkin et al. Efficient delay-tolerant particle filtering
CN102932264A (zh) 流量溢出的判断方法和装置
Shu et al. Adaptive CUSUM procedures with Markovian mean estimation
Wei et al. History-based throughput prediction with Hidden Markov Model in mobile networks
US10445444B2 (en) Flow rate prediction device, mixing ratio estimation device, method, and computer-readable recording medium
Borodina et al. Application of splitting to failure estimation in controllable degradation system
Mohammadpour et al. Selecting the best flood flow frequency model using multi-criteria group decision-making
CN104679939A (zh) 一种飞机设计经济可承受性评估过程的多准则决策方法
Wu et al. Software reliability modeling based on SVM and virtual sample
WO2023139640A1 (ja) 情報処理装置および情報処理方法
Bladt et al. Simple simulation of diffusion bridges with application to likelihood inference for diffusions
Begin et al. A DFO technique to calibrate queueing models
Sadre et al. Fitting heavy-tailed HTTP traces with the new stratified EM-algorithm
Rehman et al. Is there self-similarity in cloud QoS data?
Gunasekera Inferences on the common scale parameter of several exponential populations based on the generalized variable method

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 13877559

Country of ref document: EP

Kind code of ref document: A1

WWE Wipo information: entry into national phase

Ref document number: 14775485

Country of ref document: US

ENP Entry into the national phase

Ref document number: 2015505088

Country of ref document: JP

Kind code of ref document: A

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 13877559

Country of ref document: EP

Kind code of ref document: A1