CN103049669A - Method for forecasting water bloom and analyzing factors on basis of multivariate cyclostationary time sequence analysis and grey theory - Google Patents

Abstract

The invention discloses a method for forecasting water bloom and analyzing factors on the basis of multivariate cyclostationary time sequence analysis and a grey theory, and belongs to the technical field of environmental engineering. The method includes acquiring and preprocessing monitoring data; screening characteristic factors; modeling a time sequence of the characteristic factors; and forecasting the time sequence of the characteristic factors and analyzing the characteristic factors. A multivariate cyclostationary property of the multivariate time sequence is inspected, a process for screening the characteristic factors of the water bloom is provided, necessary conditions for reasonably modeling the multivariate time sequence of the characteristic factors of the water bloom are provided, and a basis is also provided for analyzing the characteristic factors of the water bloom. In addition, a procedure for forming the water bloom due to the multiple characteristic factors is comprehensively described, and accordingly the accuracy and the reliability of a water bloom forecasting result are improved. Besides, closeness degrees of relations among the various influencing factors and the water bloom can be obtained according to errors among the water bloom forecasting result given by the method and measured data.

Description

Wawter bloom prediction and factor approach based on the analysis of Multivariate Periodic Stationary Time Series and gray theory

Technical field

The present invention relates to a kind of wawter bloom prediction and factor approach, belong to field of environment engineering technology.

Background technology

Wawter bloom refers to appear in the eutrophication water, and when possessing suitable illumination, water temperature, weather and the hydrology etc. and be conducive to the environmental baseline of algal grown and gathering, the alga eruption sexual reproduction is assembled and also reached certain density a kind of phenomenon.Large-scale wawter bloom harm is very large, and it not only can destroy the 26S Proteasome Structure and Function of water ecosystem, and harm humans is healthy, utilizes usefulness but also can reduce water resource, threatens the sustainable exploitation and use of water resource, causes huge economic loss.In general, also lack technology and the means that effectively to administer in a short time wawter bloom at present.Therefore before wawter bloom is effectively administered, the Accurate Prediction that carries out of wawter bloom is convenient to relevant department and is taked counter-measure, harm reduction, therefore, the wawter bloom control has important scientific meaning and using value.

In the preventing and controlling of wawter bloom, the wawter bloom prediction always is a difficult point, this is because the stochastic process that wawter bloom forms is the physics of a complexity, chemistry and biological synthesis course of reaction, to form relevant characteristic factor more with wawter bloom, interaction relationship is close, consist of the aquatic ecosystem of a complexity, because the complicacy of ecological process inherent mechanism, the impact of mankind's activity and water quality information automatic telemetry, the restriction of the conditions such as reception, the origin cause of formation to wawter bloom formation mechanism and outburst is also fully unclear at present, thereby directly sets up relatively difficulty of wawter bloom mechanism ecological model.Along with going deep into of wawter bloom forecasting research, some scholars is attempted setting up accurately reliable mathematical model and is come the stochastic process that wawter bloom forms is predicted, makes some progress.Yet existing wawter bloom Forecasting Methodology still exists precision of prediction not high, the problem that prediction step is short.

In addition, the degree of correlation that different characteristic factor and wawter bloom form also is not quite similar, so the factor analysis of wawter bloom formation also is a vital task in the wawter bloom preventing and controlling.

Time series analysis is a kind of important mathematical prediction instrument of studying stochastic process.The basic assumption of tradition monobasic time series forecasting is that available information is included in the historical data, historical data has determined the variation of future time sequence, depend on the time series historical data but this hypothesis is too much, and ignored relevance between the things, also ignored the impact of external factor on the time series tendency.The model and forecast of multivariate time series has then not only been considered the impact of time series historical data on future, more consider multiple covariance and the related coefficient vital role in modeling process, having more standby application potential than monobasic time series modeling and prediction, is a kind of method that is suitable for describing and predicting wawter bloom formation stochastic process under the effect of various features factor.Therefore adopt polynary Time Series Analysis Method, the polynary time series modeling of characteristic factor to wawter bloom forms is predicted as a kind of effective way thereby carry out wawter bloom.

Yet although the abundant randomness information in the mining data of Time Series Analysis Method, the method exists the assurance of precision of prediction to depend on the limitation of large sample amount.If adopt conventional DIRECT FORECASTING METHOD to carry out forecasting of time series model, then precision of prediction is difficult to guarantee, thereby has reduced the confidence level that wawter bloom predicts the outcome.Gray theory be research minority according to probabilistic theory, namely study minority according to probabilistic background under, the decision-making of the foundation of the processing of data, the analysis of phenomenon, model, the prediction of development trend, things and the control of system and the assessment of state.Adopt the polynary temporal model of Grey Theory Forecast, can on the basis that keeps the time series modeling advantage, overcome well the limitation of forecasting of time series model, improve the precision of prediction of model.

In addition, polynary time series analysis not only is suitable for modeling and forecasting, also is a kind of effective tool of factor analysis.Therefore can utilize polynary time series analysis to predict the outcome the various characteristic factors in the wawter bloom forming process are carried out factor analysis.

Use that polynary time series analysis and gray theory carry out wawter bloom prediction and factor analysis also needs to solve following problem:

Since wawter bloom have randomness, the randomness fluctuation of each characteristic factor in the wawter bloom forming process is described, not only need to have considered the cross correlation between a plurality of characteristic factors of synchronization, to consider that also each characteristic factor is at difference autocorrelation constantly.

2. the generation of wawter bloom is relevant with seasonal variations, thereby the generation of wawter bloom also has periodically, the cyclic fluctuation of each characteristic factor in the wawter bloom forming process is described, need not only to consider that periodic environment seasonal variations on the impact of wawter bloom forming process, also will consider the reciprocal effect of multiple periodicity variation between the different characteristic factor.

3. the randomness of characteristic factor and the reciprocal effect of cyclic fluctuation between the different characteristic factor also are the problems that needs consideration in the wawter bloom forming process.

4. how will be applicable to the few data forecast Grey System Method and combine with the temporal model of building, thereby the precision of prediction of the raising temporal model of building needs specifically research.

5. how to carry out the analysis of wawter bloom characteristic factor based on the forecasting of time series model result who builds.

Summary of the invention

The objective of the invention is to predict the outcome not accurate enough in order to solve existing wawter bloom, reciprocal effect modeling difficulty in the wawter bloom forming process between the various features factor, and the decision problem of the degree of relevancy of different characteristic factor and wawter bloom generation, take the technological means based on the analysis of Multivariate Periodic Stationary Time Series and gray theory, reach wawter bloom prediction effect and factor analysis result by the characteristic factor time series modeling analysis in the wawter bloom forming process being obtained more conform to actual conditions.

For ease of explanation, all unexplained nouns and alphabetical implication are explained by following hypothesis in this instructions: the characteristic factor relevant with the wawter bloom phenomenon is divided into two kinds: a kind of is to affect the characteristic factor that wawter bloom occurs, such as nitrogen, phosphorus, pH value, dissolved oxygen DO, water temperature, illuminance etc. below is called influence factor; Another kind is to characterize the characteristic factor that wawter bloom occurs, and such as chlorophyll concentration, algae density etc. below is called the sign factor.With Y _tExpression t characteristic factor vector constantly; With y _ItRepresent i characteristic factor at t value constantly, total sampling time is N, t=1, and 2 ..., N, total n characteristic factor, i=1,2 ..., n.

Wawter bloom prediction and factor approach based on the analysis of Multivariate Periodic Stationary Time Series and gray theory provided by the invention mainly comprises following five steps:

Step 1, Monitoring Data collection and pre-service;

A plurality of characteristic factors that may affect or characterize the wawter bloom generation are monitored.Therefore the Monitoring Data temporal evolution of a plurality of characteristic factors that collected by monitoring equipment and changing is multivariate time series.Generally, primitive character factor Monitoring Data includes abnormal sample point, average drifting and part sampled point shortage of data etc. phenomenon, is difficult to directly it be carried out time series analysis.The impact that time series analysis is caused for eliminating these phenomenons, the fitting precision of raising characteristic factor temporal model, and the initial value of unified primitive character factor sequential are tackled the original sequential of each characteristic factor and are done respectively pre-service.

Step 2, characteristic factor screening;

Wawter bloom have randomness, and relevant with seasonal variations, thereby the generation of wawter bloom also has periodically, therefore the polynary sequential of its characteristic factor presents randomness, cyclical variation at occurring in nature, but the change mechanism between its a plurality of characteristic factors is metastable, can temporal evolution, i.e. the variance of the polynary sequential of characteristic factor time to time change not, so the characteristic factor sequential should be the Multivariate Periodic Stationary Time Series.

Need the polynary sequential of pretreated characteristic factor is carried out the Multivariate Periodic stationary test for this reason, and carry out characteristic factor according to Multivariate Periodic stationary test result and screen, the characteristic factor that only has the Multivariate Periodic stationary test to pass through just is fit to adopt the inventive method.The present invention is based on delay crosscorrelation matrix method and provide Multivariate Periodic stationary test method, the step 2 in the visible embodiment.Characteristic factor screening technique of the present invention is as follows:

1, the monobasic Cyclostationarity of each characteristic factor is analyzed.Because the monobasic cyclic stationary is stably necessary condition of Multivariate Periodic, therefore at first the monobasic cycle sequential of each characteristic factor is carried out the check of monobasic Cyclostationarity.Adopt the autocorrelation function method of inspection.The characteristic factor that the check of monobasic Cyclostationarity is not passed through is rejected;

All sign factors of 2, the monobasic Cyclostationarity being upchecked, such as chlorophyll concentration, algae density etc., the Multivariate Periodic stationary test method check of adopting the present invention to provide, finding out can be by the combination of comprising of Multivariate Periodic stationary test of maximum sign factors;

3, to the sign factor combination by the Multivariate Periodic stationary test, add the influence factor by the check of monobasic stability of period as much as possible, forming the combination of New Characteristics factor, until till the new characteristic factor combination that forms can not or comprise whole influence factors by the Multivariate Periodic stationary test.

Like this, according to the characteristic factor the selection result, following steps three and step 4 are carried out in the characteristic factor combination that meets the Multivariate Periodic smooth conditions.

Step 3, characteristic factor time series modeling;

1, determines the characteristic factor sequential organization;

The polynary sequential of wawter bloom characteristic factor presents randomness, cyclical variation at occurring in nature, but does not exist tendency to change, therefore with t characteristic factor vector Y constantly _tBe decomposed into periodic term C _tWith random entry R _tStack,

Y _t＝C _t+R _t (1)

$Y_t=\left(\begin{array}{c}{y}_{1t}\\ y_{2t}\\ .\\ .\\ .\\ y_{\mathrm{nt}}\end{array}\right),$ $C_t=\left(\begin{array}{c}{c}_{1t}\\ c_{2t}\\ .\\ .\\ .\\ c_{\mathrm{nt}}\end{array}\right),$ $R_t=\left(\begin{array}{c}{r}_{1t}\\ r_{2t}\\ .\\ .\\ .\\ r_{\mathrm{nt}}\end{array}\right)$

C wherein _ItBe the periodic term of i characteristic factor, r _ItBe the random entry of i characteristic factor, i=1,2 ..., n, wherein n is the number of characteristic factor.

2, set up characteristic factor sequential periodic term model;

Because cyclical variation exists interactivity impact, periodic term C between a plurality of characteristic factors _tCan't reflect that then this interactivity impact, the present invention's employing are applicable to the potential periodic regularity of mining data and reflect that mutual sex multiple latent periodic model is described between multiple periodicity if adopt the monobasic periodic model to describe, namely

Wherein C (t) is the multiple latent periodic function of multiple latent periodic model, and q is the cycle angular frequency number of diving, ω _jBe j angular frequency, a _IjBe j the amplitude that angular frequency is corresponding of i characteristic factor of multiple latent periodic model,

Be j the phase place that angular frequency is corresponding of i characteristic factor of multiple latent periodic model, i=1,2 ..., n.

3, set up characteristic factor sequential random entry model;

From Y _tDeduct C _tAfter, to random entry, i.e. Y _tSteady randomness part R _tAdopt the multivariate autoregressive model description that most widely used general, modeling is simple and be suitable for predicting, namely

$R_t=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{R}_{t-j}+E_t---\left(3\right)$

$E_t=\left(\begin{array}{c}{\mathrm{\ϵ}}_{1t}\\ {\mathrm{\ϵ}}_{2t}\\ .\\ .\\ .\\ {\mathrm{\ϵ}}_{\mathrm{nt}}\end{array}\right)$

Wherein p is the multivariate autoregressive exponent number, H _jBe n * n multivariate autoregressive matrix of coefficients, R _T-jBe the random entry of inscribing when the t-j, E _tBe n dimension white noise vector, η _IkjBe that i characteristic factor is to the multivariate autoregressive coefficient of k characteristic factor, ε _ItBe the white noise of i characteristic factor, i, k=1,2 ..., n.

4, set up characteristic factor Multivariate Periodic Stationary Time Series model.

With periodic term C _tMultiple latent periodic model formula (2) and random entry R _tMultivariate autoregressive model formation (3) combination.Because also may there be reciprocal effect in wawter bloom forming process complicated mechanism between the randomness of its characteristic factor and the cyclical variation, so random entry R _tWith periodic term C _tBetween also may have correlativity, if only with the simple addition of model of the two, then can ignore this correlativity.For this reason, in order to improve the models fitting precision, the present invention carries out interative computation with the multiple latent periodic model substitution random entry multivariate autoregressive model of periodic term, proposes the multiple latent cycle multivariate autoregressive mixture model of characteristic factor

$Y_t=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{Y}_{t-j}+C^{*}\left(t\right)+E_t---\left(4\right)$

Y wherein _T-jBe the characteristic factor vector of inscribing when the t-j, $C^{*}\left(t\right)=\left(\begin{array}{c}\underset{j=1}{\overset{q}{\mathrm{\Σ}}}{b}_{1j}\cos ({\mathrm{\ω}}_{j}t+{\mathrm{\φ}}_{1j})\\ \underset{j=1}{\overset{q}{\mathrm{\Σ}}}{b}_{2j}\cos ({\mathrm{\ω}}_{j}t+{\mathrm{\φ}}_{2j})\\ .\\ .\\ .\\ \underset{j=1}{\overset{q}{\mathrm{\Σ}}}{b}_{\mathrm{nj}}\cos ({\mathrm{\ω}}_{j}t+{\mathrm{\φ}}_{\mathrm{nj}})\end{array}\right)$ Be the multiple latent periodic function of multiple latent cycle multivariate autoregressive mixture model, b _IjBe j the amplitude that angular frequency is corresponding of i characteristic factor of multiple latent cycle multivariate autoregressive mixture model, φ _IjBe j the phase place that angular frequency is corresponding of i characteristic factor of multiple latent cycle multivariate autoregressive mixture model, i=1,2 ..., n, its computing method are seen the embodiment step 3.

Step 4, characteristic factor time series forecasting;

1, the multiple latent cycle multivariate autoregressive mixture model prediction of characteristic factor;

According to optimum prediction principle (being the minimum mean-square error forecast principle), the present invention provides the multiple latent cycle multivariate autoregressive mixture model of characteristic factor and predicts forward the l best predictor Y in step constantly at t _T+lComputing formula:

$Y_{t+l}=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{Y}_{t+l-j}+C^{*}(t+l)---\left(5\right)$

2, multiple latent cycle multivariate autoregressive mixture model gray prediction.

For the multiple latent cycle multivariate autoregressive mixture model of characteristic factor, process the sequential to be predicted that obtains after the check directly to future anticipation if use based on Monitoring Data, because the restriction of Monitoring Data sample size, its precision is often undesirable.System's gray prediction nesting in the gray theory is the system's grey forecasting model for certain structure, with GM(1, and 1) model insertion GM(1, N) model solution, to obtain the predicted value of each behavior variable.The present invention is based on system's gray prediction nesting, with each step predicted value of the multiple latent cycle multivariate autoregressive mixture model of characteristic factor last data as sequential to be predicted, and leave out first data in the sequential to be predicted, rebulid again the multiple latent cycle multivariate autoregressive mixture model of characteristic factor, make model constantly obtain revising, by that analogy, until predict till the desired step number, this nested predicted method allows sequential to be predicted be exhaled the old and inhale the new, and is usually high than direct precision of prediction.

Step 5, characteristic factor analysis;

Polynary time series analysis not only is suitable for prediction, also is a kind of effective tool of factor analysis.The present invention utilizes the polynary time series analysis of front to predict the outcome and provides the factor approach of wawter bloom characteristic factor.

1, predicated error is calculated;

Characteristic factor the selection result according to step 2, after adopting step 3 among the present invention and step 4 to set up multiple latent cycle multivariate autoregressive mixture model and carry out gray prediction to the Multivariate Periodic Stationary Time Series of the characteristic factor combination that meets the Multivariate Periodic smooth conditions, predicted data and the measured data comparison of the factor of sign can be obtained predicated error, and the computing method of predicated error are seen the step 5 in the embodiment.

2, factor analysis;

The principle of factor analysis is:

(1) characterize the less characteristic factor combination of factor predicated error, influence factor is higher with the correlativity of the factor of sign in its combination, otherwise then lower;

(2) in a characteristic factor combination, if after adding one or more influence factors, its characterize factor predicated error and add before compare significantly and diminish, the influence factor and the factor of sign significant correlation that add of explanation then, its conspicuousness criterion is seen the step 5 in the embodiment;

(3) in the combination of characteristic factor, if after adding one or more influence factors, it characterizes not compare before factor predicated error and influence factor add and significantly diminishes, and then the influence factor and the sign factor that add of explanation is substantially irrelevant.

By this factor analysis, can obtain qualitatively the degree of correlation of each influence factor in the wawter bloom characteristic factor and the factor of sign, i.e. the level of intimate that contacts of each influence factor and wawter bloom generation phenomenon.

The invention has the advantages that:

1. the present invention proposes based on the Multivariate Periodic stationary test method that postpones the crosscorrelation matrix method, solved the Multivariate Periodic stationary test problem of polynary sequential.

2. the present invention has provided the screening technique of wawter bloom characteristic factor according to the Multivariate Periodic stationary test method that proposes, and for the Rational Model of the polynary sequential of wawter bloom characteristic factor provides necessary condition, also provides the foundation for the factor analysis of wawter bloom characteristic factor.

3. the present invention adopts multiple latent periodic model to describe the cyclic fluctuation of each characteristic factor in the wawter bloom forming process, considered that periodic environment seasonal variations is on the impact of wawter bloom forming process, it is compared with the latent periodic model of monobasic, has also considered the reciprocal effect of multiple periodicity variation between the different characteristic factor.

4. the present invention adopts the randomness fluctuation of each characteristic factor in the multivariate autoregressive model description wawter bloom forming process, it is compared with the monobasic autoregressive model, not only considered the cross correlation between a plurality of characteristic factors of synchronization, considered that simultaneously each characteristic factor is at difference autocorrelation constantly.

5. the present invention proposes the multiple latent cycle multivariate autoregressive mixture model of wawter bloom characteristic factor, this model is with the interative computation of multiple latent periodic model substitution multivariate autoregressive model, and two kinds of model parameters are united find the solution, thereby randomness and the reciprocal effect of cyclic fluctuation between the different characteristic factor of characteristic factor in the wawter bloom forming process have been considered, with only the method for two kinds of simple additions of model is compared, the present invention describes more comprehensive to the wawter bloom forming process of many characteristic factors, thereby has improved the accuracy that wawter bloom predicts the outcome.

6. the present invention proposes the gray prediction method of multiple latent periodic autoregressive model, overcome the limitation that multiple latent cycle multivariate autoregressive mixture model precision of prediction depends on the large sample amount, improve the precision of prediction of model, thereby improved the confidence level that wawter bloom predicts the outcome.

7. the present invention proposes the wawter bloom characteristic factor analytical approach based on polynary time series analysis, the wawter bloom that provides according to the inventive method error with measured data that predicts the outcome can obtain the level of intimate that each influence factor and wawter bloom generation phenomenon contact.

Description of drawings

Fig. 1 the present invention is based on the wawter bloom prediction of Multivariate Periodic Stationary Time Series and gray theory and the process flow diagram of factor approach;

Fig. 2 is the process flow diagram of the multiple latent cycle multivariate autoregressive mixture model gray prediction of the present invention;

Fig. 3 is the process flow diagram of the factor analysis of characteristic factor combination of the present invention;

Fig. 4 is sequential after the original sequential of pH value among the embodiment 1 and the pre-service;

Fig. 5 is sequential after the original sequential of oxygen utilization among the embodiment 1 and the pre-service;

Fig. 6 is sequential after the original sequential of water temperature among the embodiment 1 and the pre-service;

Fig. 7 is sequential after the original sequential of turbidity among the embodiment 1 and the pre-service;

Fig. 8 is sequential after the original sequential of ammonia nitrogen among the embodiment 1 and the pre-service;

Fig. 9 is sequential after the original sequential of total nitrogen among the embodiment 1 and the pre-service;

Figure 10 is sequential after the original sequential of total phosphorus among the embodiment 1 and the pre-service;

Figure 11 is sequential after the original sequential of dissolved oxygen DO among the embodiment 1 and the pre-service;

Figure 12 is sequential after the original sequential of embodiment 1 Determination of Chlorophyll and the pre-service;

Figure 13 is sequential after the original sequential of algae density among the embodiment 1 and the pre-service;

Figure 14 adopts four kinds of methods to chlorophyllous predicted data and measured data among the embodiment 1;

Figure 15 adopts four kinds of methods to predicted data and the measured data of algae density among the embodiment 1;

The curve numbering is respectively among the figure: the original sequential of each individual features factor among the 1-embodiment 1, sequential after the pre-service of each individual features factor among the 2-embodiment 1, the measured data of each corresponding sign factor among the 3-embodiment 1, each corresponding sign factor adopts the predicted data of monobasic cyclic stationary time series analysis among the 4-embodiment 1, each corresponding sign factor adopts the predicted data of monobasic cyclic stationary time series analysis and gray theory among the 5-embodiment 1, the predicted data that each corresponding sign factor adopts the Multivariate Periodic Stationary Time Series to analyze among the 6-embodiment 1, each corresponding sign factor adopts the predicted data of the analysis of Multivariate Periodic Stationary Time Series and gray theory among the 7-embodiment 1.

Embodiment

The present invention is described in further detail below in conjunction with drawings and Examples 1.

The present invention is a kind of wawter bloom prediction and factor approach based on the analysis of Multivariate Periodic Stationary Time Series and gray theory.

For ease of explanation, all unexplained nouns and alphabetical implication are explained by following hypothesis in this instructions: the characteristic factor relevant with the wawter bloom phenomenon is divided into two kinds: a kind of is to affect the characteristic factor that wawter bloom occurs, such as nitrogen, phosphorus, pH value, dissolved oxygen DO, water temperature, illuminance etc. below is called influence factor; Another kind is to characterize the characteristic factor that wawter bloom occurs, and such as chlorophyll concentration, algae density etc. below is called the sign factor.With Y _tExpression t characteristic factor vector constantly; With y _ItRepresent i characteristic factor at t value constantly, total n characteristic factor, i=1,2 ..., n.

The concrete grammar implementing procedure is realized as shown in Figure 1 as follows:

Step 1, Monitoring Data collection and pre-service;

A plurality of characteristic factors that may affect or characterize the wawter bloom generation are monitored.Generally, the primitive character factor Monitoring Data that is collected by monitoring equipment includes abnormal sample point, average drifting and part sampled point shortage of data etc. phenomenon, tackles the original sequential of each characteristic factor and does respectively pre-service.Characteristic factor Monitoring Data preprocess method of the present invention is:

1, remove the abnormal sample point, the present invention removes absolute value greater than the abnormal sample point of 3 times of original sequential standard deviations.

2, remove average drifting, the present invention deducts its average with original sequential, and making the sequential average after the removal is zero.

3, fill up missing data, the present invention with the sampled point of the previous sampled data of the sampled point of missing data and missing data after the average of a sampled data replace the sampled point of missing data.

Step 2, characteristic factor screening;

The polynary sequential of pretreated characteristic factor is carried out the Multivariate Periodic stationary test, and carry out the characteristic factor screening according to Multivariate Periodic stationary test result.Its Multivariate Periodic stationary test method is as follows:

Remember one group of characteristic factor vector sequential Y to be tested _t, t=1,2 ..., the covariance matrix of N is

Γ ₀＝E[(Y _t-E(Y _t))(Y _t-E(Y _t)) ^T]

Γ then ₀Be n * n square formation.Y _tThe crosscorrelation matrix be

y ₀＝D ^-1Γ ₀D ^-1

Wherein

By Γ ₀Diagonal element (Γ ₁₁(0) ..., Γ _Nn(0)) forms.Y then ₀Middle element can be expressed as

$y_{\mathrm{ij}}\left(0\right)=\frac{{\mathrm{\Γ}}_{\mathrm{ij}}\left(0\right)}{\sqrt{{\mathrm{\Γ}}_{\mathrm{ii}}\left(0\right){\mathrm{\Γ}}_{\mathrm{jj}}\left(0\right)}},i,j=\mathrm{1,2},\·\·\·,n$

Γ wherein _Ij(0) expression Γ ₀In the capable j column element of i.

With Y _tThe covariance matrix on delay k rank be designated as

Γ _k＝E[(Y _t-E(Y _t))(Y _t-k-E(Y _t)) ^T],t＝1,2,…,N

Y wherein _T-kExpression t-k characteristic factor vector constantly.Γ _kBe similarly n * n square formation.Y _tDelay k rank crosscorrelation matrix be

y _k＝D ^-1Γ _kD ^-1

Y then _kMiddle element is

$y_{\mathrm{ij}}\left(k\right)=\frac{{\mathrm{\Γ}}_{\mathrm{ij}}\left(k\right)}{\sqrt{{\mathrm{\Γ}}_{\mathrm{ii}}\left(0\right){\mathrm{\Γ}}_{\mathrm{jj}}\left(0\right)}},i,j=\mathrm{1,2},\·\·\·,n$

Wherein, Γ _Ij(k) expression Γ _kIn the capable j column element of i.

Γ _kCan by

${\mathrm{\Γ}}_k=\frac{1}{N-k}\underset{t=k+1}{\overset{N}{\mathrm{\Σ}}}(Y_t-\stackrel{\‾}{Y}){(Y_{t-k}-\stackrel{\‾}{Y})}^T$

Approximate evaluation.Wherein $\stackrel{\‾}{Y}=\frac{1}{N}\underset{t=k+1}{\overset{N}{\mathrm{\Σ}}}{Y}_t.$

Work as Y _tDuring for the Multivariate Periodic sequential, establishing its cycle is T _i, i=1,2 ..., q then has

$Y_t=Y_{t+T_i}$

Y then _tDelay k+T _iThe covariance matrix on rank is

${\mathrm{\Γ}}_{k+T_i}=E\left[(Y_t-E\left(Y_t\right)){(Y_{t-(k+T_i)}-E\left(Y_t\right))}^T\right]=E\left[(Y_t-E\left(Y_t\right)){(Y_{t-k}-E\left(Y_t\right))}^T\right]={\mathrm{\Γ}}_k$

As can be known, work as Y _tDuring for the Multivariate Periodic sequential, it postpones the covariance matrix Γ on k rank _kAlso have periodically.

If so front k y=(y ₀, y ₁..., y _k), k≤N-1 has the character of cyclic fluctuation, then Y _tThe Multivariate Periodic stationary test pass through.

Step 3, characteristic factor time series modeling;

1, determines the characteristic factor sequential organization;

With t characteristic factor vector Y constantly _tBy formula (1) is decomposed into periodic term C _tWith random entry R _tStack.

2, set up characteristic factor sequential periodic term model;

The present invention adopts multiple latent periodic model formula (2) Expressive Features factor sequential periodic term C _t

3, set up characteristic factor sequential random entry model;

From Y _tDeduct C _tAfter, be Y to random entry _tSteady randomness part R _tAdopt multivariate autoregressive model formation (3) to describe.

4, set up characteristic factor Multivariate Periodic Stationary Time Series model;

With periodic term C _tMultiple latent periodic model formula (2) and random entry R _tMultivariate autoregressive model formation (3) be combined into multiple latent cycle multivariate autoregressive mixture model formula (4).Its associated methods is as follows:

If matrix coefficient polynomial expression

$H\left(z\right)=I-\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{z}^j$

Wherein z is the matrix coefficient root of polynomial, and I is n * n unit matrix.

Then according to formula (3), have

$H\left(B\right)R_t=(I-\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{B}^j)R_t=R_t-\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{R}_{t-j}=E_t$

Wherein B is the backward Shift operators of time t.

The following formula both sides are with taking advantage of H ^-1(B), random entry R then _tCan be expressed as

R _t＝H ^-1(B)Ε _t

And have

$\det (I-\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{z}^j)\≠0,\left|z\right|\≤1$

With periodic term C _tThe real-valued form formula (2) of multiple latent periodic model write as the form of complex value

$C_t=C\left(t\right)=\left(\begin{array}{c}\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{1j}\exp \left(i{\mathrm{\λ}}_{j}t\right)\\ \underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{2j}\exp \left(i{\mathrm{\λ}}_{j}t\right)\\ .\\ .\\ .\\ \underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{\mathrm{nj}}\exp \left(i{\mathrm{\λ}}_{j}t\right)\end{array}\right)=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)---\left(6\right)$

α wherein _KjBe j the amplitude that angular frequency is corresponding of k characteristic factor of multiple latent periodic model complex value form, λ _jBe j angular frequency of multiple latent periodic model complex value form, ${\mathrm{\α}}_j=\left(\begin{array}{c}{\mathrm{\α}}_{1j}\\ {\mathrm{\α}}_{2j}\\ .\\ .\\ .\\ {\mathrm{\α}}_{\mathrm{nj}}\end{array}\right).$

Utilize formula (6) formula (1) can be write as

$Y_t=C_t+R_t=\left(\begin{array}{c}\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{1j}\exp \left(i{\mathrm{\λ}}_{j}t\right)\\ \underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{2j}\exp \left(i{\mathrm{\λ}}_{j}t\right)\\ .\\ .\\ .\\ \underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{\mathrm{nj}}\exp \left(i{\mathrm{\λ}}_{j}t\right)\end{array}\right)+H^{-1}\left(B\right)E_t=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)+H^{-1}\left(B\right)E_t---\left(7\right)$

Be multiplied by simultaneously H (B) on following formula (7) both sides, get Η ₀=-I obtains

$H\left(B\right)Y_t=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}H\left(B\right){\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)+E_t$

$=-\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}\left(\underset{l=0}{\overset{p}{\mathrm{\Σ}}}{H}_l{\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_j(t-l)\right)\right)+E_t$

$=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}(-\underset{l=0}{\overset{p}{\mathrm{\Σ}}}{H}_l{\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_{j}l\right))\exp \left(i{\mathrm{\λ}}_{j}t\right)+E_t$

$=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\β}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)+E_t$

Wherein

β _Kj≠ 0, k=1,2 ..., n represents j the amplitude that angular frequency is corresponding of k characteristic factor of multiple latent cycle multivariate autoregressive mixture model complex value form.

So, formula (7) can be write as

$Y_t=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{Y}_{t-j}+\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\β}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)+E_t---\left(8\right)$

Formula (8) is the complex value form of multiple latent cycle multivariate autoregressive mixture model.The real-valued form of its correspondence is multiple latent cycle multivariate autoregressive mixture model formula (4).

The amplitude b of the multiple latent cycle multivariate autoregressive mixture model of real-valued form _IjAnd phase _Ij, i=1,2 ..., n, j=1,2 ..., q is estimated as:

If λ _q=π gets b _Iq=β _Iq, φ _Iq=0, ω _q=π,

If λ _j∈ (0, π), get b _Ij=2| β _Ij|, ω _j=λ _j, φ _Ij=arg (β _Ij).

Step 4, characteristic factor time series forecasting;

1, characteristic factor Multivariate Periodic Stationary Time Series model prediction;

According to the minimum mean-square error forecast principle, obtain the multiple latent cycle multivariate autoregressive mixture model of characteristic factor and predict forward the l minimum mean-square error forecast formula (5) in step constantly at t.

2, characteristic factor Multivariate Periodic Stationary Time Series gray prediction;

Employing system gray prediction nesting to the flow process of characteristic factor Multivariate Periodic Stationary Time Series prediction as shown in Figure 2, concrete steps are as follows:

(1) set prediction total step number H, make u=1,2 ..., H represents the step number predicted.

(2) establish { Y _Ut, t=1,2 ..., one group of characteristic factor Multivariate Periodic Stationary Time Series to be predicted when N represents u step prediction is set up multiple latent periodic model

C wherein _Ut, C _u(t), q _u, a _Uij, ω _Uj,

I=1,2 ..., n is the periodic term when carrying out u step prediction, multiple latent periodic function and the corresponding model parameter of multiple latent periodic model.

(3) from Y _UtIn deduct periodic term C _Ut, remaining data are set up the multivariate autoregressive model

$R_{\mathrm{ut}}=\underset{j=1}{\overset{p_u}{\mathrm{\Σ}}}{H}_{\mathrm{uj}}{R}_{u(t-j)}+E_{\mathrm{ut}}---\left(10\right)$

R wherein _Ut, R _{U (t-j)}, p _u, Η _UjRandom entry and corresponding model parameter when being u step prediction.Ε _UtWhite noise when being u step prediction.

(4) set up multiple latent cycle multivariate autoregressive mixture model

$Y_{\mathrm{ut}}=\underset{j=1}{\overset{p_u}{\mathrm{\Σ}}}{H}_{\mathrm{uj}}{Y}_{u(t-j)}+C_u^{*}\left(t\right)+E_{\mathrm{ut}}---\left(11\right)$

P wherein _u, Η _Uj,

Model parameter when being u step prediction and the multiple latent periodic function of multiple latent cycle multivariate autoregressive mixture model.

(5) adopt the minimum mean-square error forecast formula (5) of multiple latent cycle multivariate autoregressive mixture model to Y _UNCarry out forward one-step prediction

$Y_{u(N+1)}=\underset{j=1}{\overset{p_u}{\mathrm{\Σ}}}{H}_{\mathrm{uj}}{Y}_{u(N+1-j)}+C_u^{*}(N+1)=Y_{N+u}---\left(12\right)$

Y wherein _N+uCharacteristic factor is vectorial constantly for N+u.So far finish the one-step prediction of characteristic factor Multivariate Periodic Stationary Time Series.

(6) remove first data Y of this group sequential _U1, the one-step prediction value Y of adding sequential _{U (N+1)}, obtain one group of new sequential { Y _Ut, t=2,3 ..., N+1, the sequential { Y to be predicted when this sequential was predicted as the u+1 step _{(u+1) t}, t=1,2 ..., N.

(7) make u=u+1, then temporal representation to be predicted is { Y _Ut, t=1,2 ..., N repeats above-mentioned steps (2)～(6), can dope one by one the polynary sequential { Y of characteristic factor _N+u, u=1,2 ..., H.

Step 5, characteristic factor analysis;

1, predicated error is calculated;

Characteristic factor the selection result according to step 2, after adopting step 3 among the present invention and step 4 to set up multiple latent cycle multivariate autoregressive mixture model and carry out gray prediction to the Multivariate Periodic Stationary Time Series of the characteristic factor combination that meets the Multivariate Periodic smooth conditions, predicted data and its measured data comparison with the moment of the factor of sign can be obtained predicated error, the expression mode of predicated error is the average error absolute value, namely

Wherein, H is the total number of predicted data.

2, factor analysis;

The flow process of the factor analysis of characteristic factor combination as shown in Figure 3, concrete steps are:

(1) computational representation factor predicated error makes up the sign factor as analyzing combination.

(2) choosing arbitrarily an influence factor of not analyzing from characteristic factor combination joins and analyzes in the combination, again the sign factor of analyzing in the combination is calculated predicated error, if significantly diminish with phase ratio error before the adding, this influence factor then is described and characterizes the factor analysis height, it is retained in analyzes combination.Otherwise illustrate that this influence factor is low with the sign factor analysis, rejects it from analyzing the combination.Its conspicuousness criterion is

Add after this influence factor before predicated error decrease 〉=this influence factor of adding 1% of predicated error

If contain a plurality of sign factors in the characteristic factor combination, then carry out conspicuousness with a plurality of sign factor predicated error sums and judge.

(3) repeat above-mentioned steps (2), until equal analyzed the going over of all influence factors in the characteristic factor combination, then having analyzed the influence factor that keeps in the combination is with the sign factor or with wawter bloom related close influence factor occurs, the influence factor of analyzing in the combination is added the size ordering that changes characterizing the factor predicated error by it, obtain qualitatively the degree of correlation of each influence factor and the factor of sign.

Embodiment 1:

Step 1, characteristic factor Monitoring Data gather and pre-service;

10 wawter bloom characteristic factors to Taihu Lake, Jiangsu Province year June in June, 2009 to 2012 are monitored, and specifically see Table 1.

Table 1 wawter bloom characteristic factor monitoring list

Title

The pH value

Oxygen utilization

Water temperature

Turbidity

Ammonia nitrogen

Total nitrogen

Total phosphorus

Dissolved oxygen DO

Chlorophyll

Algae density

Unit

Nothing

mg/L

℃

NTU

mg/L

mg/L

mg/L

mg/L

mg/L

Individual/L

These two characteristic factors of its Determination of Chlorophyll and algae density are the sign factor, and remaining 8 characteristic factor is influence factor.Monitoring equipment has recorded altogether 1104 days wawter bloom characteristic factor data, and its 10 original sequential of characteristic factor are seen respectively the Grey curves 1 among Fig. 4 to Figure 13.To the pre-service that the original sequential of each characteristic factor is removed respectively the abnormal sample point, removed average drifting and fill up the sampled point of missing data, pretreated characteristic factor sequential is seen the black curve 2 among Fig. 4 to Figure 13.

Step 2, characteristic factor screening;

The polynary sequential of pretreated characteristic factor is carried out the Multivariate Periodic stationary test, and carry out the characteristic factor screening according to Multivariate Periodic stationary test result.The characteristic factor the selection result is as follows:

1, adopt the autocorrelation function method of inspection to carry out the check of monobasic Cyclostationarity to all characteristic factors, the check of the monobasic Cyclostationarity of turbidity, ammonia nitrogen, these three characteristic factors of total phosphorus is not passed through as can be known, does not consider so it should be removed.

The sign factor of 2, the monobasic Cyclostationarity being upchecked, i.e. chlorophyll concentration and algae density, the Multivariate Periodic stationary test method that adopts the present invention to provide, these the two kinds combinations that characterize factors are by check.

3, to characterizing the combination of factor, add respectively the influence factor by the check of monobasic stability of period, be pH value, oxygen utilization, water temperature, total nitrogen, these 5 influence factors of dissolved oxygen DO forming the combination of New Characteristics factor, the characteristic factor that is comprised of these 5 influence factors and 2 sign factors as can be known through check makes up and can pass through the Multivariate Periodic stationary test.

Step 3, characteristic factor time series modeling;

Get front 1095 days characteristic factor sequential and be used for modeling, thereby can predict the 1096th day to the 1104th day these data of 9 days according to front 1095 day data.

1, determines the characteristic factor sequential organization;

With the characteristic factor sequential by formula (1) be decomposed into the stack of periodic term and random entry.

2, set up characteristic factor sequential periodic term model;

Adopt multiple latent periodic model formula (2) Expressive Features factor sequential periodic term.

3, set up characteristic factor sequential random entry model;

After deducting periodic term from the characteristic factor sequential, adopt multivariate autoregressive model formation (3) to describe to random entry.

4, set up characteristic factor Multivariate Periodic Stationary Time Series model;

Be multiple latent cycle multivariate autoregressive mixture model formula (4) with the multiple latent periodic model of periodic term and the multivariate autoregressive models coupling of random entry.

Step 4, the polynary time series forecasting of wawter bloom characteristic factor;

1, characteristic factor Multivariate Periodic Stationary Time Series model prediction;

According to the minimum mean-square error forecast principle, obtain the minimum mean-square error forecast formula (5) of the multiple latent cycle multivariate autoregressive mixture model of characteristic factor.

2, characteristic factor Multivariate Periodic Stationary Time Series gray prediction;

Employing system gray prediction nesting is as follows to 1104 days 9 steps prediction concrete steps of front 1095 days characteristic factor time series forecasting to the:

(1) set prediction total step number 9, make u=1,2 ..., the step number of 9 expression predictions.

(2) establish { Y _Ut, t=1,2 ..., one group of characteristic factor Multivariate Periodic Stationary Time Series to be predicted during 1095 expression u step prediction is set up multiple latent periodic model formula (9).

(3) from Y _UtIn remove periodic term C _Ut, remaining data are set up multivariate autoregressive model formation (10).

(4) set up multiple latent cycle multivariate autoregressive mixture model formula (11).

(5) adopt the minimum mean-square error forecast formula (5) of multiple latent cycle multivariate autoregressive mixture model to Y _U1095Carry out forward one-step prediction formula (12), obtain 1095+u days characteristic factor time series forecasting value Y _1095+uSo far finish the one-step prediction of characteristic factor sequential.

(6) remove first data Y of this group sequential _U1, the one-step prediction value Y of adding sequential _U1096, obtain one group of new sequential { Y _Ut, t=2,3 ..., 1096, the sequential { Y to be predicted when this sequential was predicted as the u+1 step _{(u+1) t}, t=1,2 ..., 1095.

(7) make u=u+1, then temporal representation to be predicted is { Y _Ut, t=1,2 ..., 1095, repeat above-mentioned steps (2)～(6), can dope one by one the polynary sequential { Y of characteristic factor _1095+u, u=1,2 ..., 9, i.e. the 1096th day to the 1104th day characteristic factor data.

Step 5, characteristic factor analysis;

1, predicated error is calculated;

According to the characteristic factor the selection result of step 2, behind the characteristic factor combination employing the inventive method modeling and forecasting that meets the Multivariate Periodic smooth conditions, it characterizes factor predicted data and the measured data comparison can obtain corresponding predicated error.

2, factor analysis is as shown in table 2.

Table 2 sign factor prediction average error absolute value

As can be known from Table 2, these two predicated errors that characterize factor of chlorophyll and algae density are respectively 0.3857 and 1.3910, and the predicated error sum is 1.7767; After adding the influence factor dissolved oxygen DO, two predicated errors that characterize factor are respectively 0.3055 and 1.3153, and its predicated error sum has reduced 0.1559 before comparing adding, and decrease is greater than 1% of predicated error sum before adding, so can keep; After adding influence factor pH value, characterizing has increased by 0.0038 before factor predicated error sum is compared adding on the contrary, so should reject; After adding the influence factor water temperature, before comparing adding, sign factor predicated error sum altogether reduced 0.4398, and decrease is greater than 1% of predicated error sum before adding, so can keep; After adding the influence factor total nitrogen, before comparing adding, sign factor predicated error sum altogether reduced 0.3191, and decrease is greater than 1% of predicated error sum before adding, so can keep; After adding the influence factor oxygen utilization, before comparing adding, sign factor predicated error sum altogether reduced 0.04, and decrease is less than 1% of predicated error sum before adding, so should reject.

Therefore, the final influence factor that keeps is dissolved oxygen DO, water temperature, total nitrogen, illustrate that these two degrees of correlation that characterize factor of these three influence factors and chlorophyll and algae density are higher, the relation that is the generation of these three influence factors and wawter bloom is the closest, and these three influence factors by its predicated error decrease from big to small, and the level of intimate that also namely occurs with wawter bloom sorts from high to low and is followed successively by water temperature, total nitrogen, dissolved oxygen DO.

In addition, be the advantage that illustrates that more intuitively the present invention adopts the analysis of Multivariate Periodic Stationary Time Series and gray theory that wawter bloom is predicted, the present invention provides in addition and adopts the time series analysis of monobasic cyclic stationary, the time series analysis of monobasic cyclic stationary and gray theory, the analysis of Multivariate Periodic Stationary Time Series to predict the outcome to the polynary sequential of wawter bloom characteristic factor, to compare with the inventive method result.

Adopt four kinds of methods that the 1096th day to the 1104th day the polynary sequential of wawter bloom characteristic factor is predicted, its 2 characterize factor predicted data and measured data respectively such as Figure 14 and shown in Figure 15, and the prediction average error absolute value of its 7 characteristic factors is as shown in table 3.

Four kinds of method predictions of table 3 average error absolute value

As seen, compare other three kinds of methods based on the polynary time series forecasting result of the wawter bloom characteristic factor of the inventive method and more conform to, predict that the average error absolute value is less with measured result.

Claims (5)

1. based on wawter bloom prediction and the factor approach of the analysis of Multivariate Periodic Stationary Time Series and gray theory, may further comprise the steps:

Step 1, Monitoring Data collection and pre-service;

Described pre-service comprises removal abnormal sample point, removes average drifting and fill up the sampled point missing data;

Step 2, carry out characteristic factor screening according to Multivariate Periodic stationary test method;

Step 3, characteristic factor time series modeling;

(1), determines the characteristic factor sequential organization;

With t characteristic factor vector Y constantly _tBe decomposed into periodic term C _tWith random entry R _tStack,

Y _t＝C _t+R _t (1)

$Y_t=\left(\begin{array}{c}{y}_{1t}\\ y_{2t}\\ .\\ .\\ .\\ y_{\mathrm{nt}}\end{array}\right),$ $C_t=\left(\begin{array}{c}{c}_{1t}\\ c_{2t}\\ .\\ .\\ .\\ c_{\mathrm{nt}}\end{array}\right),$ $R_t=\left(\begin{array}{c}{r}_{1t}\\ r_{2t}\\ .\\ .\\ .\\ r_{\mathrm{nt}}\end{array}\right)$

C wherein _ItBe the periodic term of i characteristic factor, r _ItBe the random entry of i characteristic factor, i=1,2 ..., n, wherein n is the number of characteristic factor;

(2), set up characteristic factor sequential periodic term model;

Adopt multiple latent periodic model to describe, namely

Wherein C (t) is the multiple latent periodic function of multiple latent periodic model, and q is the cycle angular frequency number of diving, ω _jBe j angular frequency, a _IjBe j the amplitude that angular frequency is corresponding of i characteristic factor of multiple latent periodic model,

Be j the phase place that angular frequency is corresponding of i characteristic factor of multiple latent periodic model, i=1,2 ..., n;

(3), set up characteristic factor sequential random entry model;

From Y _tDeduct C _tAfter, to random entry, i.e. Y _tSteady randomness part R _tAdopt the multivariate autoregressive model description, namely

$R_t=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{R}_{t-j}+E_t---\left(3\right)$

$E_t=\left(\begin{array}{c}{\mathrm{\ϵ}}_{1t}\\ {\mathrm{\ϵ}}_{2t}\\ .\\ .\\ .\\ {\mathrm{\ϵ}}_{\mathrm{nt}}\end{array}\right)$

Wherein p is the multivariate autoregressive exponent number, H _jBe n * n multivariate autoregressive matrix of coefficients, R _T-jBe the random entry of inscribing when the t-j, E _tBe n dimension white noise vector, η _IkjBe that i characteristic factor is to the multivariate autoregressive coefficient of k characteristic factor, ε _ItBe the white noise of i characteristic factor, i, k=1,2 ..., n;

(4), set up characteristic factor Multivariate Periodic Stationary Time Series model:

The multiple latent periodic model substitution random entry multivariate autoregressive model of periodic term is carried out interative computation, propose the multiple latent cycle multivariate autoregressive mixture model of characteristic factor

$Y_t=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{Y}_{t-j}+C^{*}\left(t\right)+E_t---\left(4\right)$

Y wherein _T-jBe the characteristic factor vector of inscribing when the t-j, $C^{*}\left(t\right)=\left(\begin{array}{c}\underset{j=1}{\overset{q}{\mathrm{\Σ}}}{b}_{1j}\cos ({\mathrm{\ω}}_{j}t+{\mathrm{\φ}}_{1j})\\ \underset{j=1}{\overset{q}{\mathrm{\Σ}}}{b}_{2j}\cos ({\mathrm{\ω}}_{j}t+{\mathrm{\φ}}_{2j})\\ .\\ .\\ .\\ \underset{j=1}{\overset{q}{\mathrm{\Σ}}}{b}_{\mathrm{nj}}\cos ({\mathrm{\ω}}_{j}t+{\mathrm{\φ}}_{\mathrm{nj}})\end{array}\right)$ Be the multiple latent periodic function of multiple latent cycle multivariate autoregressive mixture model, b _IjBe j the amplitude that angular frequency is corresponding of i characteristic factor of multiple latent cycle multivariate autoregressive mixture model, φ _IjBe j the phase place that angular frequency is corresponding of i characteristic factor of multiple latent cycle multivariate autoregressive mixture model, i=1,2 ..., n;

Step 4, characteristic factor time series forecasting;

(1), the multiple latent cycle multivariate autoregressive mixture model prediction of characteristic factor;

According to the minimum mean-square error forecast principle, provide the multiple latent cycle multivariate autoregressive mixture model of characteristic factor and predict forward the l best predictor Y in step constantly at t _T+lComputing formula:

$Y_{t+l}=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{Y}_{t+l-j}+C^{*}(t+l)---\left(5\right)$

(2), multiple latent cycle multivariate autoregressive mixture model gray prediction;

Based on system's gray prediction nesting, with each step predicted value of the multiple latent cycle multivariate autoregressive mixture model of characteristic factor last data as sequential to be predicted, and leave out first data in the sequential to be predicted, rebulid again the multiple latent cycle multivariate autoregressive mixture model of characteristic factor, make model constantly obtain revising, by that analogy, until predict till the desired step number;

Step 5, characteristic factor analysis obtain the degree of correlation of each influence factor in the wawter bloom characteristic factor and the factor of sign qualitatively, i.e. the level of intimate that contacts of each influence factor and wawter bloom generation phenomenon.

2. the wawter bloom based on the analysis of Multivariate Periodic Stationary Time Series and gray theory according to claim 1 predicts and factor approach that the stationary test of Multivariate Periodic described in step 2 method is as follows:

Remember one group of characteristic factor vector sequential Y to be tested _t, t=1,2 ..., the covariance matrix of N is

Γ ₀＝E[(Y _t-E(Y _t))(Y _t-E(Y _t)) ^T]

Γ then ₀Be n * n square formation, Y _tThe crosscorrelation matrix be

y ₀＝D ^-1Γ ₀D ^-1

Wherein

By Γ ₀Diagonal element (Γ ₁₁(0) ..., Γ _Nn(0)) forms, then y ₀Middle element representation is

$y_{\mathrm{ij}}\left(0\right)=\frac{{\mathrm{\Γ}}_{\mathrm{ij}}\left(0\right)}{\sqrt{{\mathrm{\Γ}}_{\mathrm{ii}}\left(0\right){\mathrm{\Γ}}_{\mathrm{jj}}\left(0\right)}},i,j=\mathrm{1,2},\·\·\·,n$

Γ wherein _Ij(0) expression Γ ₀In the capable j column element of i;

With Y _tThe covariance matrix on delay k rank be designated as

Γ _k＝E[(Y _t-E(Y _t))(Y _t-k-E(Y _t)) ^T],t＝1,2,…,N

Y wherein _T-kExpression t-k characteristic factor vector constantly, Γ _kBe similarly n * n square formation.Y _tDelay k rank crosscorrelation matrix be

y _k＝D ^-1Γ _kD ^-1

Y then _kMiddle element is

$y_{\mathrm{ij}}\left(k\right)=\frac{{\mathrm{\Γ}}_{\mathrm{ij}}\left(k\right)}{\sqrt{{\mathrm{\Γ}}_{\mathrm{ii}}\left(0\right){\mathrm{\Γ}}_{\mathrm{jj}}\left(0\right)}},i,j=\mathrm{1,2},\·\·\·,n$

Wherein, Γ _Ij(k) expression Γ _kIn the capable j column element of i;

Γ _kBy

${\mathrm{\Γ}}_k=\frac{1}{N-k}\underset{t=k+1}{\overset{N}{\mathrm{\Σ}}}(Y_t-\stackrel{\‾}{Y}){(Y_{t-k}-\stackrel{\‾}{Y})}^T$

Approximate evaluation, wherein $\stackrel{\‾}{Y}=\frac{1}{N}\underset{t=k+1}{\overset{N}{\mathrm{\Σ}}}{Y}_t;$

Work as Y _tDuring for the Multivariate Periodic sequential, establishing its cycle is T _i, i=1,2 ..., q then has

$Y_t=Y_{t+T_i}$

Y then _tDelay k+T _iThe covariance matrix on rank is

${\mathrm{\Γ}}_{k+T_i}=E\left[(Y_t-E\left(Y_t\right)){(Y_{t-(k+T_i)}-E\left(Y_t\right))}^T\right]=E\left[(Y_t-E\left(Y_t\right)){(Y_{t-k}-E\left(Y_t\right))}^T\right]={\mathrm{\Γ}}_k$

Work as Y _tDuring for the Multivariate Periodic sequential, it postpones the covariance matrix Γ on k rank _kAlso have periodically;

If so front k y=(y ₀, y ₁..., y _k), k≤N-1 has the character of cyclic fluctuation, then Y _tThe Multivariate Periodic stationary test pass through.

3. the wawter bloom based on the analysis of Multivariate Periodic Stationary Time Series and gray theory according to claim 1 predicts and factor approach that the associated methods of the multiple latent cycle multivariate autoregressive mixture model of characteristic factor is as follows in the step 3:

If matrix coefficient polynomial expression

$H\left(z\right)=I-\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{z}^j$

Wherein z is the matrix coefficient root of polynomial, and I is n * n unit matrix;

Then according to formula (3), have

$H\left(B\right)R_t=(I-\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{B}^j)R_t=R_t-\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{R}_{t-j}=E_t$

Wherein B is the backward Shift operators of time t;

The following formula both sides are with taking advantage of H ^-1(B), random entry R then _tBe expressed as

R _t＝H ^-1(B)Ε _t

And have

$\det (I-\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{z}^j)\≠0,\left|z\right|\≤1$

With periodic term C _tThe real-valued form formula (2) of multiple latent periodic model write as the form of complex value

$C_t=C\left(t\right)=\left(\begin{array}{c}\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{1j}\exp \left(i{\mathrm{\λ}}_{j}t\right)\\ \underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{2j}\exp \left(i{\mathrm{\λ}}_{j}t\right)\\ .\\ .\\ .\\ \underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{\mathrm{nj}}\exp \left(i{\mathrm{\λ}}_{j}t\right)\end{array}\right)=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)---\left(6\right)$

α wherein _KjBe j the amplitude that angular frequency is corresponding of k characteristic factor of multiple latent periodic model complex value form, λ _jBe j angular frequency of multiple latent periodic model complex value form, ${\mathrm{\α}}_j=\left(\begin{array}{c}{\mathrm{\α}}_{1j}\\ {\mathrm{\α}}_{2j}\\ .\\ .\\ .\\ {\mathrm{\α}}_{\mathrm{nj}}\end{array}\right);$

Utilize formula (6) that formula (1) is write as

$Y_t=C_t+R_t=\left(\begin{array}{c}\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{1j}\exp \left(i{\mathrm{\λ}}_{j}t\right)\\ \underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{2j}\exp \left(i{\mathrm{\λ}}_{j}t\right)\\ .\\ .\\ .\\ \underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_{\mathrm{nj}}\exp \left(i{\mathrm{\λ}}_{j}t\right)\end{array}\right)+H^{-1}\left(B\right)E_t=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)+H^{-1}\left(B\right)E_t---\left(7\right)$

Be multiplied by simultaneously H (B) on following formula (7) both sides, get Η ₀=-I obtains

$H\left(B\right)Y_t=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}H\left(B\right){\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)+E_t$

$=-\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}\left(\underset{l=0}{\overset{p}{\mathrm{\Σ}}}{H}_l{\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_j(t-l)\right)\right)+E_t$

$=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}(-\underset{l=0}{\overset{p}{\mathrm{\Σ}}}{H}_l{\mathrm{\α}}_j\exp \left(i{\mathrm{\λ}}_{j}l\right))\exp \left(i{\mathrm{\λ}}_{j}t\right)+E_t$

$=\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\β}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)+E_t$

Wherein

β _Kj≠ 0, k=1,2 ..., n represents j the amplitude that angular frequency is corresponding of k characteristic factor of multiple latent cycle multivariate autoregressive mixture model complex value form;

So, formula (7) is write as

$Y_t=\underset{j=1}{\overset{p}{\mathrm{\Σ}}}{H}_j{Y}_{t-j}+\underset{j=1}{\overset{2q}{\mathrm{\Σ}}}{\mathrm{\β}}_j\exp \left(i{\mathrm{\λ}}_{j}t\right)+E_t---\left(8\right)$

Formula (8) is the complex value form of multiple latent cycle multivariate autoregressive mixture model, and the real-valued form of its correspondence is multiple latent cycle multivariate autoregressive mixture model formula (4);

The amplitude b of the multiple latent cycle multivariate autoregressive mixture model of real-valued form _IjAnd phase _Ij, i=1,2 ..., n, j=1,2 ..., q is estimated as:

If λ _q=π gets b _Iq=β _Iq, φ _Iq=0, ω _q=π,

If λ _j∈ (0, π), get b _Ij=2| β _Ij|, ω _j=λ _j, φ _Ij=arg (β _Ij).

4. the wawter bloom based on the analysis of Multivariate Periodic Stationary Time Series and gray theory according to claim 1 is predicted and factor approach, adopts system's gray prediction nesting as follows to the concrete steps of characteristic factor Multivariate Periodic Stationary Time Series prediction in the step 4:

(1) set prediction total step number H, make u=1,2 ..., H represents the step number predicted;

(2) establish { Y _Ut, t=1,2 ..., one group of characteristic factor Multivariate Periodic Stationary Time Series to be predicted when N represents u step prediction is set up multiple latent periodic model

C wherein _Ut, C _u(t), q _u, a _Uij, ω _Uj, I=1,2 ..., n is the periodic term when carrying out u step prediction, multiple latent periodic function and the corresponding model parameter of multiple latent periodic model;

(3) from Y _UtIn deduct periodic term C _Ut, remaining data are set up the multivariate autoregressive model

$R_{\mathrm{ut}}=\underset{j=1}{\overset{p_u}{\mathrm{\Σ}}}{H}_{\mathrm{uj}}{R}_{u(t-j)}+E_{\mathrm{ut}}---\left(10\right)$

R wherein _Ut, R _{U (t-j)}, p _u, Η _UjRandom entry and corresponding model parameter when being u step prediction.Ε _UtWhite noise when being u step prediction;

(4) set up multiple latent cycle multivariate autoregressive mixture model

$Y_{\mathrm{ut}}=\underset{j=1}{\overset{p_u}{\mathrm{\Σ}}}{H}_{\mathrm{uj}}{Y}_{u(t-j)}+C_u^{*}\left(t\right)+E_{\mathrm{ut}}---\left(11\right)$

P wherein _u, Η _Uj,

Model parameter when being u step prediction and the multiple latent periodic function of multiple latent cycle multivariate autoregressive mixture model;

(5) adopt the minimum mean-square error forecast formula (5) of multiple latent cycle multivariate autoregressive mixture model to Y _UNCarry out forward one-step prediction

$Y_{u(N+1)}=\underset{j=1}{\overset{p_u}{\mathrm{\Σ}}}{H}_{\mathrm{uj}}{Y}_{u(N+1-j)}+C_u^{*}(N+1)=Y_{N+u}---\left(12\right)$

Y wherein _N+uCharacteristic factor is vectorial constantly for N+u; So far finish the one-step prediction of characteristic factor Multivariate Periodic Stationary Time Series;

(6) remove first data Y of this group sequential _U1, the one-step prediction value Y of adding sequential _{U (N+1)}, obtain one group of new sequential { Y _Ut, t=2,3 ..., N+1, the sequential { X to be predicted when this sequential was predicted as the u+1 step _{(u+1) t}, t=1,2 ..., N;

(7) make u=u+1, then temporal representation to be predicted is { Y _Ut, t=1,2 ..., N repeats above-mentioned steps (2)～(6), dopes one by one the polynary sequential { Y of characteristic factor _N+u, u=1,2 ..., H.

5. the wawter bloom based on the analysis of Multivariate Periodic Stationary Time Series and gray theory according to claim 1 predicts and factor approach that the concrete steps of characteristic factor analysis are in the step 5:

(1) the sign factor predicated error in the combination of calculated characteristics factor makes up the sign factor as analyzing combination;

(2) choosing arbitrarily an influence factor of not analyzing from characteristic factor combination joins and analyzes in the combination, again the sign factor of analyzing in the combination is calculated predicated error, if significantly diminish with phase ratio error before the adding, this influence factor then is described and characterizes the factor analysis height, it is retained in analyzes combination; Otherwise illustrate that this influence factor is low with the sign factor analysis, rejects it from analyzing the combination; Its conspicuousness criterion is:

Add after this influence factor before predicated error decrease 〉=this influence factor of adding 1% of predicated error

If contain a plurality of sign factors in the characteristic factor combination, then carry out conspicuousness with a plurality of sign factor predicated error sums and judge;

(3) repeat above-mentioned steps (2), until equal analyzed the going over of all influence factors in the characteristic factor combination, then having analyzed the influence factor that keeps in the combination is with the sign factor or with wawter bloom related close influence factor occurs, the influence factor of analyzing in the combination is added the size ordering that changes characterizing the factor predicated error by it, obtain qualitatively the degree of correlation of each influence factor and the factor of sign.

Images