WO2022239609A1 - モジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラム - Google Patents

モジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラム Download PDF

Info

Publication number
WO2022239609A1
WO2022239609A1 PCT/JP2022/018118 JP2022018118W WO2022239609A1 WO 2022239609 A1 WO2022239609 A1 WO 2022239609A1 JP 2022018118 W JP2022018118 W JP 2022018118W WO 2022239609 A1 WO2022239609 A1 WO 2022239609A1
Authority
WO
WIPO (PCT)
Prior art keywords
prediction
pairwise
model
output
input
Prior art date
Application number
PCT/JP2022/018118
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 株式会社 東芝
Publication of WO2022239609A1 publication Critical patent/WO2022239609A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G05—CONTROLLING; REGULATING
    • G05B—CONTROL OR REGULATING SYSTEMS IN GENERAL; FUNCTIONAL ELEMENTS OF SUCH SYSTEMS; MONITORING OR TESTING ARRANGEMENTS FOR SUCH SYSTEMS OR ELEMENTS
    • G05B23/00—Testing or monitoring of control systems or parts thereof
    • G05B23/02—Electric testing or monitoring

Definitions

  • Infrastructure systems such as water and sewage systems, rainwater drainage systems, power systems, and transportation systems, or process-based industrial plants such as steel processes, petrochemical processes, semiconductor manufacturing processes, etc., usually measure multiple process conditions.
  • Multiple online sensors are installed.
  • a system called a process monitoring and control system SCADA: Supervisory Control And Data Acquisition
  • SCADA Supervisory Control And Data Acquisition
  • process data flow rate, temperature, water quality, operation amount, etc.
  • Process monitoring and control systems usually provide these time-series data as trend graphs to observers in many cases.
  • Examples of infrastructure systems in which plant monitoring using prediction information is used include water demand prediction in water purification systems, sewage inflow prediction to sewage treatment plants, and rainwater inflow to sewage treatment plants and rainwater drainage pumping stations during rainy weather. Forecasting, forecasting river water levels and dam water levels in rainy weather, short-term rainfall forecasts (nowcast) using weather radar, etc., traffic volume forecasts and congestion forecasts in transportation systems, power demand forecasts in power systems, buildings, etc. There are systems such as predicting the number of customers in cafeterias and stores in facilities, and they are widely used in various fields in the infrastructure field.
  • a plant monitoring system that uses forecast information not only provides the forecast information itself as added value, but also uses forecast information for plant operation planning and plant control, resulting in more efficient operation and safer operation. It can also be used to achieve high performance. For example, demand forecasts such as water demand and electricity demand can be used to achieve efficient operations by optimizing water operation plans and power generation and storage plans (EMS: Energy Management System). can.
  • EMS Energy Management System
  • rainfall prediction and inflow prediction contribute to risk reduction such as avoidance and suppression of flooding by using prediction values as input information for operation control of rainwater drainage pumps.
  • a modular time-series data prediction device will be described as a basic means for performing adjustable time-series data prediction.
  • the modular time-series prediction device of this embodiment allows a prediction model (system) to be decomposed into partial models (subsystems), and each partial model can be easily attached and detached.
  • the modular time-series data prediction device for example, data for any one of multiple monitoring variables is missing, or multiple monitoring variables include outliers (abnormal data). If the monitored variable is unstable measurement data containing a large amount of noise or the like, it is possible to prevent the influence of the data from spreading to the entire system.
  • the modular time-series data prediction device for example, it is possible to appropriately determine whether or not explanatory variables are necessary (unnecessary), easily detach a partial model, or partially adjust the prediction model. You can
  • the advantage is that the model is "decomposable" due to its modular configuration.
  • the advantage of adopting a modular configuration in the prediction model is not only the advantage obtained simply by being decomposable. It makes it possible to build a convincing model, which in turn facilitates adjustments.
  • An identifiable problem is a problem when the structure of the model (the form of the formula expressing the model) is limited in advance, for example, when the model is limited to a linear transfer function model or a time series model, the parameters included in the model ) can be uniquely determined. In this case, it is particularly called a parameter identifiability problem among the identifiability problems.
  • This parameter identifiability problem is a problem that frequently appears even in problems that deal with real data, and the following three cases are particularly frequent problems.
  • the first is the problem of multicollinearity mentioned above, which is a problem when there is a strong correlation between input variables (explanatory variables). For example, when predicting rainwater inflow to a rainwater pumping station, radar rainfall information is used as an input variable, and the rainfall amounts of two different meshes that are close to each other are used as input variables. It is a problem that appears in , and in actual data, it is extremely common that there is a strong correlation between explanatory variables.
  • PE condition Persistent Exciting Condition
  • y(t) b1 ⁇ u(t-1)+b2 ⁇ u(t-2)
  • u(t-1) u
  • (t-2) k (constant value)
  • b1 and b2 cannot be determined uniquely.
  • a generalization of this concept is a parameter identifiable condition called a continuous excitation condition, but such a condition is often not satisfied.
  • the above three examples are typical cases in which the identifiability of parameters is lacking, but in actual problems, such situations are encountered very often.
  • the object of modeling is a relatively small mechanical system or electrical system, it is common practice to conduct experiments to improve identifiability and build a model after securing identifiability.
  • the data to be modeled is data obtained from natural phenomena (rainfall, solar radiation, wind speed, etc.), or the target plant
  • the data is in the operating state, and it is not possible to conduct experiments to improve the identifiability of the parameters, and the data that can be used for modeling is given. Because there are so many, it is often inherently difficult to avoid the lack of identifiability.
  • the model must be built without identifiability.
  • the relationship between u1 and y and the relationship between u2 and y are determined by simple regression using equations (3) and (4) to construct two prediction models. Then, the average of the outputs y of the two prediction models is calculated. In this way, following a constructive procedure, the values of a1 and a2 are each estimated as 0.25, which is a physically plausible and convincing value, and as a result, (1 ), the coefficients of the multiple regression equation can be defined.
  • XAI Explainable AI
  • the so-called XAI is a top-down approach in which the whole system is treated as a single large system and a black box model is constructed, and then the model is decomposed into appropriately interpretable forms and explained later. (or an inductive approach).
  • the concept of the modular time-series data prediction device of this embodiment is to integrate models (subsystems) that can be rationally explained that can be understood relatively easily individually. It is a bottom-up (or deductive) approach to building predictive models.
  • the above bottom-up approach cannot consider mutual interference between subsystems in advance, so when focusing on prediction accuracy, it may be inferior to the top-down (XAI) approach. cannot be denied.
  • the bottom-up approach has the potential to be superior to the top-down approach in terms of rational explainability, as it enables model construction that maintains the identifiability of parameters as explained earlier. high.
  • White-box modeling is essentially a deductive, constructive modeling technique.
  • a model is constructively constructed using physically meaningful parameters. is directly related to the input/output to identify parameters that are easy to interpret, and construct a model using them.
  • the modular time-series data prediction device of this embodiment is based on the above insight.
  • FIG. 1 is a diagram schematically showing a rainfall inflow prediction system to which a modular time-series data prediction device of one embodiment is applied.
  • the target process 1 of the rainfall inflow prediction system is a flowmeter 11, a trunk flowmeter 12, a trunk water level gauge 131-13K, a ground rain gauge 141-14M, a radar rain gauge 15, and pump well water level gauges 161, 162. , an inflow culvert 17 , a storm pump well 18 , a storm pump 19 and an inflow gate 110 .
  • the radar rain gauge 15 includes radar rain gauges 1511-15QP for each mesh of the Q ⁇ P mesh.
  • the various sensors 11-15, 161, and 162 of the urban rainwater drainage process 1 measure quantities representing the state of the process and quantities related to driving operations at predetermined intervals (eg, 30 seconds, 60 seconds).
  • the values measured by various sensors 11-15, 161, 162 are collected and stored in a modular time-series data prediction device.
  • the inflow gate 110 is provided in the channel between the inflow culvert 17 and the rainwater pump well 18, and operates to adjust the amount of rainwater flowing into the rainwater pump well 18 from the inflow culvert 17 by opening and closing.
  • Water level gauges 161 and 162 are installed before and after the inflow gate 110, and the operation of the inflow gate 110 is controlled according to the values measured by the water gauges 161 and 162, for example.
  • FIG. 2 is a diagram schematically showing one configuration example of a modular time-series data prediction device according to one embodiment.
  • the modular time-series data prediction device of this embodiment includes a data collection and storage unit 2, a data extraction unit 3, a pairwise prediction model identification unit 4, a prediction model synthesis method definition unit 5, and a pairwise output variable prediction unit 6. , a synthetic output variable prediction unit 7 , a prediction error evaluation unit 8 , a pairwise prediction model correction unit 9 , and an output prediction result observation unit 10 .
  • the configuration of the modular time-series data prediction device of the present embodiment is, for example, an arithmetic device including at least one processor such as a CPU or MPU and a memory in which a program executed by the processor is stored. can be realized by software or by a combination of software and hardware.
  • the data collection and storage unit 2 collects and stores the values of process variables measured by various sensors 11-15, 161, and 162.
  • the data collection and storage section 2 includes an output variable data selection section 21 and an input variable data selection section 22 .
  • the output variable data selector 21 selects at least one of the values measured by the various sensors 11-15, 161, 162 as an output variable.
  • the input variable data selector 22 selects a plurality of input variables from the values measured by the various sensors 11-15, 161 and 162. FIG. Note that the number of input variables may be less than or equal to the number of process variables.
  • the data extraction unit 3 extracts predetermined data from the time-series data stored in the data collection storage unit 2.
  • the data extraction unit 3 includes an offline prediction model identification data extraction unit 31 , an online prediction data extraction unit 32 , and an evaluation data extraction unit 33 .
  • the data extraction unit 31 for offline prediction model identification uses data for prediction model identification from the data collection and storage unit 2 at a predetermined identification cycle or at an external request. Data for a predetermined period (before the current time) is extracted.
  • the online prediction data extraction unit 32 extracts prediction data necessary for online prediction from the measurement values of the various process sensors 11-15, 161, and 162 from the data collection and storage unit 2 at a predetermined prediction cycle in real time. Extract.
  • the evaluation data extraction unit 33 extracts the measured values of the various process sensors 11-15, 161, and 162 required for prediction error evaluation from the data collection/storage unit 2 in real time at a predetermined evaluation interval or upon request from the outside. Extract.
  • the pairwise prediction model identification unit 4 selects the output variables selected by the output variable data selection unit 21 and the input variable data selection unit 22 using the measured value data extracted by the offline prediction model identification data extraction unit 31
  • the obtained input variables are pairwise for each variable, and a prediction model is identified by the method described later.
  • the prediction model synthesis method definition unit 5 is a method of synthesizing an output variable (each output variable if there are multiple output variables) from the output value (pairwise prediction value) of the pairwise prediction model identified by the pairwise prediction model identification unit 4.
  • the synthesized output variable prediction unit 7 outputs the predicted values of the output variables synthesized according to the method defined by the prediction model synthesis method definition unit 5 for the pairwise predicted values output from the pairwise output variable prediction unit 6 .
  • a prediction error evaluation unit 8 stores a prediction value obtained by synthesizing the pairwise prediction value output from the pairwise output variable prediction unit 6 and the pairwise prediction value output from the combined output variable prediction unit 7, and evaluates the prediction error at a predetermined period. Alternatively, at a predetermined timing, the errors between the actual values of the output variables extracted by the evaluation data extracting unit 33 and the stored pairwise and combined predicted values are evaluated respectively.
  • the pairwise prediction model correction unit 9 determines whether or not adjustment of the pairwise prediction model is necessary, and corrects the necessary pairwise prediction model.
  • the pairwise prediction model correction unit 9 can, for example, adjust the parameters of the pairwise prediction models, delete the pairwise prediction models, or integrate the pairwise prediction models.
  • the pairwise prediction model correction unit 9 supplies the corrected pairwise prediction model to the pairwise output variable prediction unit 6 .
  • the output prediction result observation unit 10 includes display means such as a monitor and a user interface, and can present at least the prediction values of the synthetic output variables by the synthetic output variable prediction unit 7 to the observer.
  • the output prediction result observation unit 10 may be, for example, a personal computer, a tablet terminal, or a portable terminal such as a smart phone. If necessary, the output prediction result observation unit 10 may further present to the observer at least one of a plurality of outputs of the pairwise prediction model by the pairwise output variable prediction unit 6 and prediction distribution information synthesized therefrom. .
  • the output variable data selector 21 selects a variable to be predicted from process variables measured by various process sensors 11-15, 161, and 162, and sets the selected variable as an output variable. For example, in the urban rainwater drainage process 1, when the inflow into the normal rainwater pump well 18 is to be predicted, the output variable data selection unit 21 is measured by the flow meter (rainwater pump well inflow meter) 11. The rainwater pump well inflow is selected as the output variable to be predicted.
  • the output variable data selection unit 21 calculates and stores the calculated value of the inflow, and selects the calculated value of the inflow as the output variable.
  • the output variable data selection unit 21 may select multiple items as output variables. Also in the urban rainwater drainage process 1, for example, when rainwater flows into the rainwater pump well 18 from a plurality of locations, the output variable data selection unit 21 selects the inflow amount at each of the plurality of inflow locations as an output variable. You may
  • the input variable data selection unit 22 selects measured values that may affect the output variables selected by the output variable data selection unit 21 as input variable candidates. If there is no prior information as to whether or not there is an influence, the input variable data selector 22 may mechanically set all process variables as input variables.
  • the input variable data selection unit 22 may select variables selected as output variables as input variables. For example, in time-series data analysis, so-called autoregression, which predicts future values from past values of variables themselves selected as output variables, is often considered. In such a case, the rainwater pump well inflow by the flow meter (rainwater pump well inflow meter) 11 may be selected as both an output variable and an input variable.
  • the input variable data selection unit 22 may select another output variable as an input variable for a certain output variable. By doing so, it becomes possible to consider the correlation between the plurality of output variables.
  • the input variable data selector 22 can exclude another output variable as an input variable for a given output variable.
  • the input variable data selection unit 22 also selects other variables, for example, variables that are clearly unrelated to the output variable from information such as the physical connection relationship of pipes, and causal relationships (input (cause) and output (result)). ) and ) can be excluded from the input variables in advance.
  • the rainwater pump well water level measured by the rainwater pump well water level gauges 161 and 162 changes under the influence of the rainwater pump well inflow from the flow meter (rainwater pump well inflow meter) 11 selected as an output variable. is clear from a physical point of view, so the water level of the rainwater pump well is excluded from the input variables.
  • the input variable data selection unit 22 preliminarily performs a correlation analysis with the output variables selected in advance by the output variable data selection unit 21 for variables whose physical structural relationships are not clear, and the absolute value of the correlation coefficient is determined. Variables less than the value of (for example, 0.5) may be preliminarily excluded from the input variables. Furthermore, when performing the correlation analysis, the input variable data selection unit 22 performs the correlation analysis while shifting the time when the variable to be analyzed is measured, obtains the shifted time L when the correlation is the highest, From the sign information of L, variables having signs that do not satisfy the causal relationship may be excluded from the input variables in advance. It is highly probable that the rainwater pump well water level variable will be excluded from the input variables even by such an analysis.
  • the input variable data selector 22 uses the flow meter (rainwater pump well inflow meter) 11 excluding only the rainwater pump well water level measured by the rainwater level gauges 161 and 162.
  • the inflow of rainwater pump wells, the trunk inflow measured by the trunk flow meter 12, the trunk water level measured by the trunk water level gauges 131-13K, the ground rainfall measured by the ground rain gauges 141-14M, and Q Radar rainfall in each of the meshes 1511-15QP of the ⁇ P mesh are selected as input variables.
  • the K trunk water level gauges 131 to 13K have recently been A relatively large number of trunk water level gauges with a K value of about 10 to 100 are often installed, and the Ministry of Land, Infrastructure, Transport and Tourism has installed a weather radar called XRAIN nationwide, so there are many locations.
  • XRAIN weather radar
  • the mesh size of the Q ⁇ P mesh it is not uncommon for the mesh size of the Q ⁇ P mesh to be about 100, so the number of items in the time-series data actually measured ( number of variables) can be quite large.
  • the number of items of measurement variables is generally several hundred to several thousand, not limited to the urban rainwater drainage process 1 to which the modular time-series data prediction device of this embodiment is applied.
  • the input variables are, for example, rainwater pump well inflow by the flow meter (rainwater pump well inflow meter) 11 excluding only the measured values by the rainwater pump well water level gauges 161 and 162.
  • Measured amount measured value of trunk inflow by trunk flow meter 12, measured value of K trunk water levels by trunk water level gauges 131 to 13K, and M ground rainfall measurements by ground rain gauges 141 to 14M and the measured value of the radar rainfall of each mesh of the Q ⁇ P meshes 1511 to 15QP.
  • the number of input variables is p (positive integer).
  • the offline prediction model identification data extraction unit 31 first extracts the output variable data selection unit 21 and the input variable data selection unit 22 at the timing requested using the data for the predetermined period specified by the person who builds the prediction model. Time-series data for the corresponding predetermined period is acquired.
  • the offline prediction model identification data extraction unit 31 periodically builds a prediction model, and when it is desired to update the model periodically, the time of the input and output variables for a predetermined period with a predetermined period TL in the direction of progress of time Series data may be extracted periodically.
  • uk a column vector, and each row corresponds to each time sample of the time-series data.
  • the off-line prediction model identification data extraction unit 31 which is set as a time-series data set y of rainwater inflow, which is an output variable, has the effect of extracting the data set uk and the data set y at a predetermined cycle.
  • the pairwise prediction model identification unit 4 uses uk and y to identify prediction models for each of p pairs of (u1, y), (u2, y), ... (up, y).
  • the output variable y is only one of the rainwater pump well inflow by the flow meter (rainwater pump well inflow meter) 11, but there are a plurality (h) If there are outputs y, the same process should be performed for each of h y's. Even if the number of output variables increases to l, the subsequent process is repeated h times and the same process is performed. Without loss of generality, let us assume that there is only one output variable.
  • the most typical prediction model performs identification by assuming a model in the form of linear regression expressed by the following equation.
  • y(t) a1 ⁇ uk(tL)+a2 ⁇ uk(tL-1)+...+an ⁇ uk(tL-n+1)+c (5)
  • c is a bias parameter representing the difference between the averages of the input uk and the output y.
  • the first advantage of pairwise identification is that the concept of "multiple inputs" disappears, so the multicollinearity problem (resulting from relationships between multiple input variables) never arises. , this problem can essentially be avoided, and from the standpoint of identifiability of parameters, it is more advantageous than a multi-input single-output system. Supplement this a little.
  • the value of the input variable is not something that can be controlled (adjusted) by the constructor of the prediction model, but given observation (measurement) information. Whether or not a frequency component is included is given as an external condition. The smaller the number of parameters, the better the identifiability. Conversely, the greater the number of parameters, the worse the identifiability. more likely. Strictly speaking, the possibility of poor identifiability is 0% when the number of parameters is small. Same as identity.
  • the 1-input 1-output model of the above equation (5) includes n parameters, whereas the p-input 1-output prediction model includes pxn parameters, so the number of inputs It can be seen that the more p increases, the more likely it becomes difficult to uniquely identify the parameter.
  • the problem of multicollinearity arises, and identifiability is lost regardless of the number of parameters to be identified, and input information is Therefore, it is theoretically impossible to uniquely determine the parameter values using any method.
  • the model of formula (5) with one input and one output the problem of multicollinearity does not occur in principle. Only when frequency information is not included, in this case, the parameter ak can be made identifiable (uniquely identifiable) by adjusting the number of parameters to be identified by adjusting the order n.
  • n 1 and one parameter a1
  • the identifiability of parameters can be maintained by identifying only
  • the upper limit nmax of the number of identifiable parameters n can be obtained by examining the above sustained excitation conditions using a given uk. (Akaike Information Criterion), BIC (Bayesian Information Criterion), MDL Criterion (Minimum Coding Criterion), etc.
  • a value can be determined to always maintain the identifiability of equation (5).
  • the second advantage of pairwise identification is that even if there is only one input variable other than the output variable, this method uses an autoregressive AR model with the output variable as input and the input variable (not the output variable). This is the effect of separating and identifying the two models, the FIR model and the input variable).
  • an autoregressive component AR component
  • FIR component/MA component weighted moving average component
  • the output variable is the pump well inflow
  • the effect of the autoregressive component on the inflow is The effect of the FIR component on the ground rainfall is stronger than that of the FIR component, and even if the ground rainfall value changes suddenly due to a change in rainfall, the change in the predicted value is small until a change in the inflow appears, and the predicted value lags behind the actual measured value. phenomenon is likely to occur.
  • the pairwise method which individually identifies an autoregressive model (AR model) for inflow and a finite impulse response model (FIR model) for ground rainfall
  • AR model autoregressive model
  • FIR model finite impulse response model
  • users such as observers can determine the AR component for forecasting.
  • the influence and the influence of the FIR component can be adjusted.
  • a third advantage of pairwise identification which is closely related to the second advantage, is that the meaning of the model in equation (1) for each input variable is clear.
  • y(t) a1 ⁇ y(t-1)+a2 ⁇ y(t-2)+...+any ⁇ y(t-ny) +b1 ⁇ u(tL)+b2 ⁇ u(tL-1)+...+bnu ⁇ u(tL-nu+1) (6)
  • Equation (7) and (8) are each in the form of equation (5), so when pairwise identification is performed, Y and U each mean the predicted value of y. Therefore, it can be interpreted that Y and U are predicted values of output y for different explanatory variables (input variables).
  • pairwise identification is that a model that is easy to assign meaning can be constructed.
  • the above is the advantage of identifying a prediction model by pairwise 1-input 1-output rather than multi-input 1-output.
  • a pairwise prediction model by identifying in the form of a regression equation in the form of FIR or AR in equation (5), frequency analysis, stability analysis, etc., which are conventionally well known in the field of digital signal processing and system identification Since various analysis methods can be directly applied to pairwise, if you have knowledge in the field, you can also perform more detailed analysis and adjustment in pairwise.
  • the AR model / FIR model like the formula (5) is further decomposed, and the variables considering the time delay are decomposed into single regression models each regarded as an explanatory variable to construct a prediction model.
  • y(t) a ⁇ uk(tL)+c (9)
  • the parameters included in the equation are the delay time L, the magnification (gain, proportionality coefficient) a of the output variable with respect to the input variable uk, the input variable and the output variable
  • the parameters included in the equation are the delay time L, the magnification (gain, proportionality coefficient) a of the output variable with respect to the input variable uk, the input variable and the output variable
  • Equation (5) it is not clear and it is possible to interpret Equation (5) itself using knowledge of system identification and digital signal processing, it is difficult to intuitively interpret individual ak.
  • the interpretation of a is easy, and simply indicates the magnification of the output variable with respect to the input variable when the average value is set to 0.
  • equation (9) can be used in cases where the parameter identification results when some algorithm is applied deviate from the intuitive feeling (for example, , such things can happen when abnormal data are mixed in), for example, it becomes easier for the predictive model builder to intuitively and manually adjust the parameters, and the contents of the predictive model can be changed It will be very easy to explain to others.
  • c in equation (9) is It means the amount of inflow when it is zero, and it can be interpreted that this corresponds to the amount of sewage in the case of combined sewage where sewage and rainwater pass through the same sewer.
  • a in the formula (9) indicates the ratio of the rainwater inflow to the rainfall by the ground rain gauge when the sewage amount is subtracted from the inflow. If the correlation between the ground rainfall and the inflow is sufficiently high, , such a simple scale factor also has meaning.
  • the delay time L in equation (9) means the delay time of the flow of rainwater from the ground rainfall observation point to the rainwater pump well. It can be interpreted to mean the sum of runoff times to reach the stormwater pump well.
  • 3 and 4 are scatter diagrams schematically showing an example of the relationship between main water level data and inflow data.
  • the horizontal axis represents main water level data at a point corresponding to one of the main water level gauges 131 to 13K, and the vertical axis represents inflow data to be predicted corresponding to inflow gauge 11.
  • An example of a diagram is shown.
  • MSE mean squared error
  • Fig. 3 for example, when the value of the main water level is -1.5, the value of the inflow is predicted to be approximately 9 according to the regression line. However, according to the scatter diagram of FIG. 3, it can be seen that the inflow amount when the value of the main water level is -1.5 is overwhelmingly greater than 9.
  • the reason why there is a discrepancy between the predicted value and the measured value by the regression line is that there is an overwhelming amount of data around -2.5 to -2, so the data around this area is more suitable. This is because the coefficients are estimated as follows. Although this regression line is theoretically correct, it is not necessarily preferable in actual inflow prediction. This is because the forecasted inflow is used to control rainwater drainage pumps and support operations, so there is an implicit demand to predict as accurately as possible when the inflow increases in order to avoid flooding risks. This is because the prediction accuracy is often not so important when the inflow is small.
  • FIG. 4 shows an example of regression lines L1-L11 obtained as a result of applying the least squares method and 10 other various regression methods other than the least squares method to the same data as in FIG.
  • methods 1 to 10 are various advanced regression methods incorporated in a software package called Python, but the regression lines L1-L10 by any method are However, it was not possible to obtain an appropriate prediction value for
  • the eleventh method is a method that extracts only data in the range of main water levels from -2.2 to -1.5 and applies the ordinary least squares method without using other data. This is because the inventors of the present application understand that the regression line is obtained so as to fit the data around the main water level -2.5 to -2, and intentionally discard the data around that area. This is a method intentionally adjusted to fit the data when there is a large amount of inflow.
  • the originally desired parameters may not be obtained unless various measures are taken, such as limiting the data range or thinning out the data.
  • the straight line L11 obtained by the eleventh method is based on the extremely simple principle that a straight line can be uniquely determined by determining two points, it takes almost no effort for a person to visually draw a line like the straight line L11. don't need it.
  • manual adjustment can be performed very easily by using simple regression parameters that are extremely easy to understand and can be clearly visualized on the scatter plots shown in FIGS. That is, using a simple regression model can be interpreted using easily understood parameters of scale factor (gain) and bias, making adjustments easier.
  • the regression coefficient (proportional coefficient, slope, magnification, gain) a of the simple regression model has the following relationship with the correlation coefficient r.
  • a r ⁇ ( ⁇ y/ ⁇ u) (10)
  • r, ⁇ y, and ⁇ u are the correlation coefficient between the input variable and the output variable, the standard deviation of the output variable, and the standard deviation of the input variable, respectively.
  • the regression coefficient a and the bias c have the following relationship when identified by the method of least squares so that the prediction squared error of the simple regression model is minimized.
  • c ⁇ y-a ⁇ u (11)
  • ⁇ y and ⁇ u represent the mean of the output variables and the mean of the input variables, respectively.
  • This relationship means only the gain and bias of a simple regression model to explain the model, the mean and standard deviation of each variable (input and output variables), and the correlation coefficient between input and output.
  • the simple regression input-output relationship can be expressed using only the correlation coefficient as shown in the following equation. I can explain everything.
  • equations (9) and (10) to (12) consist only of parameters that are easy to understand, but compared to equation (5), they have an extremely simple structure. Therefore, in actual application, the expression capability is not sufficient, and there is a possibility that sufficient prediction accuracy cannot be obtained. Therefore, we consider synthesizing the form of formula (5) from the form of formula (9) (or formulas (10) to (12)).
  • a model in the form of an FIR or AR model that is formally the same as formula (1) can be obtained.
  • averaging can be simply performed as described above, but weighted averaging can be performed as follows using the relationship between the gain, bias, and correlation coefficient of the simple regression described above. can also
  • nonlinear regression such as formula (16)
  • a regression model is constructed for each input variable with an output variable, so the relationship between an input variable and an output variable is non-linear, but the relationship between another input variable and an output variable is linear.
  • (5) and (16) can be used for each input variable, and it is generally expected that this will improve the explainability.
  • the above series of actions are actions of the pairwise prediction model identification unit 4 in this embodiment.
  • the predictive model synthesizing method defining unit 5 defines a predictive model synthesizing method defined by the pairwise predictive model identifying unit 4 .
  • the simplest prediction model synthesis method is a pairwise prediction model that predicts an output variable for each input variable.
  • p prediction outputs of the pairwise prediction model identification unit 4 (hereinafter, prediction outputs for each input variable u1, u2, ..., up are y1, y2, ..., yp ) is defined as follows.
  • Equation (17) The definition of equation (17) is the most basic definition method, but by improving this definition, the reliability of the synthesized prediction output can be improved, or the prediction can be intentionally biased. It will be possible to make predictions This will be explained in order below.
  • outliers outliers, abnormal values
  • prediction accuracy not only deteriorates, but also the synthesized prediction itself may become meaningless (break down) due to being dragged by outliers.
  • the predictive model synthesis method definition unit 5 can take the following multiple approaches.
  • the first is an approach that focuses mainly on the fact that the synthesized prediction value itself does not collapse against the latter outlier, and is a method of adopting robust estimation instead of the sample averaging process of equation (17). .
  • sample mean used in equation (17) is one method of estimating the representative value from p predicted values of y1, y2, . . . , yp. In the field of statistics, this is called location parameter estimation.
  • location parameter estimation In other words, the “sample mean” is one of the methods for estimating the position parameters. It is also known to be the weakest (non-robust) in terms of robustness. Therefore, in order to improve robustness against outliers, it is better to replace Eq. (17) with robust location parameter estimation.
  • the breakdown point This is an indicator of what percentage of outliers are allowed to be mixed in the data used for statistical estimation.
  • the sample mean value becomes ⁇ just by replacing the data of one point with ⁇ , so the breakdown point of "sample mean” is 0%. is.
  • the maximum breakdown point is known to be 50% in robust statistics, and the 'median' is known as the location parameter estimator with the maximum breakdown point. Therefore, if the primary purpose is to prevent the predicted value from being affected by outliers, replace the "average" in equation (17) with the "median” and It can also be defined as
  • trimmed average is a method of averaging after deleting ⁇ % of the top and bottom of p data.
  • trimmed average is a method of averaging after deleting ⁇ % of the top and bottom of p data.
  • increasing ⁇ increases robustness
  • the limit is the median estimate
  • decreasing ⁇ to 0 decreases robustness
  • the limit is the normal average (sample mean). It is clear that Therefore, this trimmed average can be defined by the predictive model synthesis method definition unit 5 .
  • the robustness can be adjusted by adjusting ⁇ , but it may be difficult to adjust properly, so another robust estimation method can be used.
  • HL estimation Hodges-Leemann estimator
  • MCD is a method generally used for multivariate data
  • MCD is a method of extracting a predetermined percentage of p data (usually 50% to 75%) and estimating using the data that minimizes the variance among all combinations. can also be used to estimate the mean.
  • MCD is also known as a method with high robustness and good estimation efficiency.
  • HL estimation and MCD estimation require calculation for all combinations of the number of data extracted from p data (2 for HL estimation, p/2 to 3p/4 for MCD). Therefore, when the number p of input variables of the pairwise model increases, the processing time increases dramatically. there is a possibility.
  • the bootstrap method may be used instead of performing processing for "all combinations" so that processing is repeated as many times as possible in real time. That is, predetermined k pieces of data are randomly and repeatedly extracted from p pieces of y1, y2, . Among the N average values obtained repeatedly, for example, the median value of the N average values is adopted in the same manner as the HL estimation.
  • the predictive model synthesis method defining unit 5 may define a method in which robust estimation such as M-estimation, in which outliers are weighted and estimated, is applied to position parameter estimation.
  • the pairwise prediction model identification unit 4 can use the (multiple) correlation coefficient between the prediction output and the actual output when the pairwise prediction model is identified, or its squared coefficient of determination. .
  • the pairwise prediction model identification unit 4 can use the (multiple) correlation coefficient between the prediction output and the actual output when the pairwise prediction model is identified, or its squared coefficient of determination. .
  • the operation of determining the weight of the weighted average value according to the prediction ability of the pairwise prediction model and defining the prediction model synthesis method by equation (18) is the prediction model synthesis method definition unit in this embodiment. 5 is an example of the operation.
  • the third is to define the combined prediction output by the weighted average value of each pairwise prediction model in equation (17) as in the second, but this weight is identified (estimated) using data. It is determined by At this time, by using the data used in identifying the pairwise prediction model and the identified parameter values, time-series data of predicted values y1, y2, .
  • a prediction model can be synthesized by three stages of learning (estimation), that is, a transfer function model of .
  • the fourth approach is an approach when you want to have a tendency (bias) in the synthesized forecast, and the purpose is to intentionally overestimate or underestimate the forecast. This is a synthesis method of prediction output when
  • this prediction information is often used as support information for determining the timing of starting and stopping the rainwater drainage pump.
  • inflow prediction in the control of a rainwater drainage pump it is preferable to predict the excessive side in order to determine the start timing of the pump, but in order to determine the stop timing, it is preferable to predict the excessive side. Since it is preferable to use
  • the method for synthesizing prediction output based on such motivation is shown below.
  • a parameter estimation method is also defined at the same time.
  • the position parameter is the average
  • the standard deviation is adopted as the scale parameter.
  • MAD Median Absolute Deviation
  • a method of estimating the scale parameter for each position parameter is defined. In this way, the method of estimating the position parameters and scale parameters to be applied to y1, y2, . . . , yp is defined in advance. These are expressed as ⁇ and ⁇ respectively below (usually ⁇ means mean value and ⁇ means standard deviation, but here they are not necessarily the mean and standard deviation, but the position parameters and scale parameter).
  • the method of extracting K% quantiles corresponding to K% in y1, y2, . it can be written in the form of ⁇ k ⁇ , so it is also possible to specify the quantiles set in this way and the method of directly extracting them. For example, if you want to make an overestimate prediction, you can define a prediction output that is synthesized so that the prediction value corresponding to the 90% quantile is adopted and overprediction is performed, assuming that abnormal values in the prediction are excluded to some extent. It is also possible to define a synthesized prediction output that employs prediction values corresponding to the 10% quantile and performs under-prediction.
  • the prediction model synthesis method definition unit 5 can be defined as a prediction model synthesis method by any one of these or a combination thereof. Note that the series of actions described above are performed off-line, periodically or irregularly, using past data.
  • the data of the input variables are extracted from the online prediction data extraction unit 32 at the cycle TH.
  • the input variables are the rainwater pump well inflow amount by the flow meter (rainwater pump well inflow meter) 11, the trunk inflow amount by the trunk flow meter 12, and the K trunk water levels by the trunk water level gauges 131 to 13K. , M ground rainfalls from ground rain gauges 141-14M, and radar rainfalls at each of the Q ⁇ P meshes 1511-15QP.
  • the pairwise output variable prediction unit 6 uses the pairwise prediction models of each input variable and output variable identified by the pairwise prediction model identification unit 4 to calculate p prediction outputs. This is inputting the online input variable data extracted by the online prediction data extraction unit 32 to the pairwise prediction model in which the parameter values identified in the form of formulas (5), (14), (15), etc. can be calculated immediately by This is the operation of the pairwise output variable prediction section 6.
  • FIG. 1 the pairwise prediction model identification unit 4
  • the synthetic output variable prediction unit 7 using the synthesis method definition formula corresponding to the expression (17) or (18) defined in the prediction model synthesis method definition unit 5, the output from the pairwise output variable prediction unit 6 The combined prediction output is computed by inputting the pairwise prediction output for each input variable.
  • the prediction error evaluation unit 8 acquires the prediction result of the prediction output calculated by the combined output variable prediction unit 7 and the pairwise prediction output result of the pairwise output variable prediction unit 6, Save as time-series data with period TH.
  • the saved prediction result is held until the timing of error evaluation.
  • the prediction error evaluation timing may be instructed from the outside by the user of the prediction model, or may be performed at a predetermined cycle TM (for example, TH ⁇ TM ⁇ TL) in the advancing direction of time.
  • TM for example, TH ⁇ TM ⁇ TL
  • the evaluation data extraction unit 33 outputs the time series data of the prediction output stored in the prediction error evaluation unit 8 and the time-measured output variable That is, in this embodiment, time-series data of rainwater inflow is extracted.
  • the prediction error evaluation unit 8 evaluates the error between the time-series data of the extracted output (rainwater inflow) and the prediction output of the p pairwise prediction models, and the error of the combined prediction output.
  • the prediction error evaluation unit 8 based on the set error evaluation criteria, at a predetermined timing or period TM, the error between the time series data and the prediction output of the p pairwise prediction models, and the synthesized Evaluate the error with the predicted output.
  • the pairwise prediction model correction unit 9 adjusts the parameters of the pairwise prediction model, deletes a specific pairwise model, or integrates a plurality of pairwise models based on the prediction error evaluated by the prediction error evaluation unit 8. I do.
  • the pairwise prediction model correction unit 9 determines that the prediction accuracy is significantly poor (the prediction error is larger than the predetermined allowable error) based on the prediction error evaluated by the prediction error evaluation unit 8. Input-output pairs are extracted, and pairwise prediction models for the extracted input-output pairs are re-identified. Judgment of poor prediction accuracy is performed by (1) directly specifying a threshold value (allowable error) for prediction error and determining a pairwise prediction model that exceeds the threshold value as a model (or input/output pair) with poor prediction accuracy.
  • a threshold value allowable error
  • the pairwise prediction model correction unit 9 uses the most recent data of a predetermined period to re-parameter make an identification.
  • the pairwise prediction model correction unit 9 reidentifies only the parameters of the pairwise forecast model between ground rainfall and rainwater inflow.
  • the point that readjustment can be performed by partial identification in this way is a great advantage of the modular type configuration (structural approach).
  • the deletion of the second specific pairwise model can be performed by a procedure similar to the adjustment of the first parameter.
  • the pairwise prediction model correction unit 9 performs error evaluation by a method similar to the method described in the adjustment of the first parameter, and extracts a pairwise prediction model with a significantly large error and degraded accuracy. and delete the extracted pairwise model.
  • the error evaluation method is the same as the first method.
  • the integration of the third pairwise prediction model is carried out as follows.
  • the pairwise prediction model correction unit 9 compares the prediction errors of the pairwise prediction models with each other and evaluates the prediction errors of the pairwise prediction models. Predictive models are extracted and grouped as mutually similar predictive models.
  • the identifiability of parameters is basically guaranteed, so even if the prediction results are similar and the parameter values are similar, For example, the input variables should have approximately the same value, but just in case, the pairwise prediction model correction unit 9 obtains the correlation coefficient between the input variables of the pairwise prediction model extracted as a candidate, If it is equal to or greater than a threshold value (for example, a correlation coefficient of 0.95), it is determined that the pairwise prediction model should be finally integrated. Then, the pairwise prediction model correction unit 9 defines a variable obtained by averaging input variables of a plurality of pairwise prediction models determined to be finally integrated as a new combined input variable.
  • a threshold value for example, a correlation coefficient of 0.95
  • the pairwise prediction model correction unit 9 (1) adopts the average value of the previously compared coefficient vectors, or (2) re-identifies the newly synthesized averaged input variables and output variables. Determine the coefficients of the pairwise predictive model for the synthesized input variables using either method.
  • the above integration is considered to occur in the following cases, for example, in the modular time-series data prediction device of this embodiment.
  • the size of the radar rainfall mesh is approximately 250 m ⁇ 250 m in the case of the X-band weather radar called XRAIN, for example. Therefore, the rainfall data of adjacent meshes often take almost the same value.
  • the rainfall data of mutually adjacent meshes are likely to be the target of integration, and the Q ⁇ P input variables may be integrated as a small number of mutually independent synthetic input variables. is high.
  • the mesh data of radar rainfall are integrated according to the similarity of their impact on rainwater inflow, thus not only reducing the number of input variables but also improving the interpretability. things can be expected.
  • One mesh rainfall data is selected as the representative input variable, and the mesh rainfall data that is very similar to the mesh rainfall corresponding to the selected input variable will be ignored. Since the rainfall data are very similar, there is usually no problem in selecting the rainfall of the representative mesh as the input variable in this way. However, if we see localized torrential rains that have occurred frequently in recent years, for example, the meshes that happen to be selected do not have much rainfall, but the nearby meshes that were not selected as input variables have very large rainfalls. If a certain mesh is selected as the representative mesh in the event of heavy rainfall, there is a possibility of erroneously predicting that the inflow will not increase.
  • the integration of the input variables described here averages the rainfall of the meshes with similar influences, which means the average rainfall of the integrated mesh, and thus a reasonable prediction in such cases.
  • rational explainability is also improved.
  • the ability to collectively interpret the input variables to be integrated in this way and make them easier to explain is also a great advantage of the modular configuration (structural approach).
  • the predicted output calculated (synthesized) by the synthetic output variable prediction unit 7 is sent to the output prediction result observation unit 10, and is used as time-series data called a trend graph, for example, by plant managers and operators who are users. presented to.
  • the current values of real-time monitoring are displayed simultaneously with the rainwater inflow amount, which is the actual output value, and its predicted value, and the future predicted value is displayed only with the predicted output. is preferred.
  • the position parameter ⁇ y and scale parameter ⁇ y are calculated, and ⁇ y ⁇ K ⁇ y (K is set parameter), or the maximum value from the prediction output of the pairwise prediction model (or the prediction value closest to the 95% point, which is close to the maximum value considering robustness) ) or the minimum value (or the predicted value closest to the 5% point close to the minimum value considering robustness) may be displayed simultaneously as the range of the direct prediction output.
  • a mechanism is provided that enables rational interpretation that does not contradict the laws of physics for its internal parameters.
  • the basic idea for obtaining such an effect is to extract elements that may be related to the output in a situation where only the input/output data is given and the internal structure is completely unknown. , the relationship between each element and the output is defined as the internal structure of the final input/output relationship. be.
  • the modular time-series data prediction device of this embodiment can provide the same explainability and ease of adjustment as a white-box model without entering into the internal structure that cannot be known originally. is realized.
  • this modular time-series forecasting model is well known in the field of statistical machine learning. It is similar to the method called ensemble learning such as bagging and boosting described above.
  • Bagging or boosting is a prediction model called a strong learner (strong learning model) with high accuracy by collecting a lot of models called weak learners (weak learning models) that can be easily configured but not necessarily high accuracy.
  • Bagging is a method of performing strong learning on a weak learner by processing such as majority voting
  • boosting is a method of training a strong learner using the results of the weak learner.
  • many weak learners are generated by changing the data used for learning (identification). , differ in that they identify pairwise prediction models corresponding to weak learners.
  • FIG. 5 is a diagram for explaining a first example of the pairwise prediction model in the modular time series prediction device of one embodiment.
  • FIG. 5 shows an example in which the FIR model and AR model of formula (1) are used as pairwise prediction models.
  • the pairwise predictive model MA of this example comprises regression blocks A1-AN and a prediction output synthesis block A5.
  • each of the regression blocks A1-A3 outputs a weighted moving average component (FIR component/MA component) as a predicted value for the input variable.
  • the regression block AN outputs an autoregressive component (AR component) as a predicted value for an input variable (past output variable).
  • FIR component/MA component weighted moving average component
  • AR component autoregressive component
  • FIG. 6 is a diagram for explaining a second example of the pairwise prediction model in the modular time-series prediction device of one embodiment.
  • FIG. 6 shows a prediction model for intentionally over-predicting and under-predicting as an example of a prediction model used when it is desired to impart a tendency (bias) to the synthesized prediction.
  • FIG. 7 and 8 are diagrams for explaining the effect of the pairwise prediction model of the second embodiment.
  • FIG. 7 shows an example of rainwater inflow prediction values and measured values by Lasso regression with L1 regularization, which is a kind of regression method.
  • FIG. 8 shows an example of rainwater inflow predicted values and measured values by the prediction model of the second embodiment.
  • the prediction error is small, but the estimated value is slightly smaller than the measured value, indicating that the estimated value is underestimated.
  • the prediction error is small, and although the prediction accuracy is good near the peak, the predicted value at the timing when the inflow increases is lower than the measured value.
  • the predicted values can be adjusted by such overestimation and underestimation, as a result of long-term operation, the tendency (bias) that normal predicted values are underestimated or overestimated with respect to actual values , it is possible to correct the bias and perform more accurate prediction by replacing the normal predicted value with an overestimated value or an underestimated value.
  • FIG. 9 is a diagram for explaining a third example of the pairwise prediction model in the modular time series prediction device of one embodiment.
  • the AR model / FIR model such as Equation (5) is further decomposed, and the variables considering the time delay are each considered as one explanatory variable. Building predictive models.
  • the pairwise prediction model MC of this embodiment includes prediction output blocks C1-CN including a plurality of time lag blocks, a plurality of simple regression blocks, and a correlation weighted average block for each input variable 1-N, and the It further includes a prediction output value combining block C(N+1) that outputs a final combined prediction value obtained by combining the prediction values output from the plurality of prediction output blocks C1-CN.
  • the prediction output block C1 includes a plurality (m) of time delay blocks 111-11m, a plurality (m) of simple regression blocks 121-12m, and a correlation weighted average block 13.
  • a plurality of time delay blocks 111-11m are variables (eg, at times t-L11, . value).
  • a plurality of simple regression blocks 121-12m use the simple regression coefficients to calculate and output predicted values based on the variables output from the time delay blocks 111-11m.
  • the correlation-weighted average block calculates and outputs the correlation-weighted average value of the multiple predicted values output from the multiple simple regression blocks 121-12m.
  • the predictive output blocks C2-CN have the same configuration as the predictive output block C1, and therefore individual descriptions thereof will be omitted.
  • the pairwise prediction model of this embodiment can be explained by the three terms of delay, gain, and bias, and fine adjustment of parameters is easy.
  • the entire prediction model can be constructed simply by obtaining the correlation coefficient between each input variable and each output variable, making it easy to interpret the prediction model. is.
  • the pairwise prediction model of the present embodiment can obtain the same effect as the pairwise prediction model of the first embodiment.
  • FIGS. 10 and 11 are diagrams for explaining an example of the effect of the pairwise prediction model of the third embodiment.
  • FIGS. 10 and 11 show an example in which rainwater inflow prediction is performed by inputting 10 input variables into the pairwise prediction model of the third embodiment.
  • an example of the prediction result in the simplest case considering only one delay time for each input variable is shown.
  • FIG. 11 shows an example of rainwater inflow prediction results using the above estimation results. According to the prediction result of FIG. 11, it can be seen that the prediction is generally performed with a certain degree of accuracy, but the prediction delay (phase delay) is large in the portion where the inflow increases rapidly.
  • FIG. 10 shows an example of the rainwater inflow prediction result when the delay time is adjusted to 15 minutes. According to the prediction results shown in FIG. 10, it can be seen that by adjusting the delay time, the phase delay at the time of rapid increase in the inflow can be greatly reduced without significantly affecting the overall prediction accuracy.
  • the pairwise prediction model of the present embodiment it is possible to construct a prediction model only with easily interpretable parameters, and not only is it possible to improve the explainability, but it is also possible to easily fine-tune the prediction model. can.
  • the simple regression coefficient and the correlation coefficient match if the input and output data are normalized by the mean and standard deviation, respectively. can be converted to the bias and Therefore, by calculating only the average and standard deviation of each input/output data and the input/output correlation, it is possible to identify all prediction model parameters simply by specifying the delay time.
  • the calculations required to construct a prediction model are the estimation (designation) of the lag time, the average of the input/output variables, and the Only standard deviation and correlation between input and output variables, and all parameters included in the prediction model are explained using only these four types of estimated values (lag time, mean, standard deviation, correlation coefficient). becomes possible.
  • the correlation between the input variable and the output variable changes or deteriorates, or the time-series data of the input variable contains a large amount of outliers. If the accuracy is not sufficient, it is easy to partially adjust the pairwise prediction model that affects the synthesized prediction output value or remove the part that has an adverse effect without re-identifying the entire prediction model. At the same time, this makes it possible to rationally explain the reason for the degradation of prediction accuracy.
  • a modular time-series data prediction device a modular time-series data prediction method, and a program are provided that rationally explain the prediction results and can easily adjust the prediction results. be able to.

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Automation & Control Theory (AREA)
  • Testing And Monitoring For Control Systems (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

予測結果を合理的に説明するとともに、容易に予測結果を調整可能とするモジュラー型時系列データ予測装置を提供する。　実施形態による装置は、複数のプロセス変数の時系列データを所定の周期で収集保存し、複数の前記プロセス変数の中から予測対象となる出力変数を選択し、複数のプロセス変数の中から入力変数の候補を選択する機能２と、時系列データから抽出された出力変数と複数の入力変数との同定用データを用いて、１入力１出力のペアワイズ予測モデルを定義する機能４と、複数のペアワイズ予測モデルの予測値の合成法を定義する機能５と、複数のペアワイズ予測モデルに複数の入力変数の予測用データを入力して、複数の入力変数のそれぞれに対応する予測値を演算する機能６と、定義された合成法により複数の予測値を合成することにより出力変数の予測値を演算する機能７と、を有する。

Description

モジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラム
　本発明の実施形態は、モジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラムに関する。
　上下水道システム、雨水排水システム、電力システム、交通システムなどのインフラシステム、あるいは、鉄鋼プロセス、石油化学プロセス、半導体製造プロセス、などのプロセス系の産業プラントなどでは、通常、複数のプロセス状態を測定する複数のオンラインセンサが設置されている。プロセス監視制御システム（ＳＣＡＤＡ: Supervisory Control And Data Acquisition）と呼ばれるシステムは、上記のインフラシステムや産業プラントに設置されたセンサ群の計測により得られるプロセスデータ（流量、温度、水質、操作量など）を取得し、時系列データとしてサーバ上に保持している。プロセス監視制御システムは、通常、これらの時系列データをトレンドグラフとして監視員に提供することが多い。
　監視員は、上記のような直接的な時系列データのトレンド監視に加えて、予測情報を利用するプラント監視を行うことがある。予測情報を利用するプラント監視システムによれば、監視中のリアルタイムの実測値（計測値）だけでなく、監視対象の将来の値を予測して、予測結果を運転管理者に提供することにより、適切な運転計画の実現を支援することができる。予測情報を利用するプラント監視システムは、プラントの異常兆候検出および異常診断と並んで、アドバンスト監視の典型的な手段である。
　予測情報を利用するプラント監視が用いられるインフラシステムの例として、例えば、浄水システムにおける水需要の予測、下水処理場への汚水流入予測、雨天時の下水処理場および雨水排水ポンプ場への雨水流入予測、雨天時の河川水位やダム水位の予測、あるいは、気象レーダなどによる降雨の短期予測（ナウキャスト）、交通システムにおける交通量予測や渋滞予測、電力システムにおける電力需要の予測、あるいは、ビルなどの施設内の食堂や店舗などの顧客数予測、などのシステムがあり、インフラ領域の様々な分野で広く使われている。
　予測情報を利用するプラント監視システムは、予測情報そのものを付加価値として提供するだけでなく、プラントの運転計画やプラントの制御に予測情報を利用することで、より効率の良い運用やより安全性の高い運用を実現するためにも用いられ得る。例えば、水需要や電力需要等の需要予測は、水運用計画や発電および蓄電計画（ＥＭＳ：エネルギーマネージメントシステム）などの最適化に利用することで、効率の良い運用を実現するために用いることができる。また、例えば降雨予測や流入予測は、雨水排水ポンプの運転制御のための入力情報として予測値を用いることで、浸水の回避や抑制などのリスク低減に貢献する。
　インフラシステムだけでなくプロセス系のプラントである鉄鋼プロセスや石油化学プロセスなどにおいても、製品品質や歩留まりの予測などに時系列データ解析を用いることも多く、その情報は生産効率（歩留まり）向上のために利用されることも多い。
　時系列データの予測は、インフラシステムやプロセス系の産業プラントで極めて重要な役割を果たし、プラントの監視データを利用したデータドリブンの（ブラックボックス的アプローチの）予測方法として、従来から様々な方法が適用されている。
　従来からプラント監視データなどの時系列データを用いた予測技術は広く用いられているが、近年の統計的機械学習分野を中心とするＡＩ関連の理論、技術、手法の急速な進歩に伴い、時系列データ解析に対するＡＩ手法も同時並行的に進展し、分野および領域を問わず、様々なアドバンストなＡＩ手法を時系列データの予測にも適用する動きが加速している。特に、ＮＮの発展形である深層学習ネットワーク（ＤＮＮ：ディープラーニングネットワーク）などは、昨今注目を集めている代表的な方法であり、ＤＮＮはＡＩの代名詞としても使われる場合もある。
　ＡＩの代表的な手法であるＤＮＮの他にも、重回帰分析で問題になる多重共線性による不安定化（悪条件問題）を回避する方法として（先に述べたＰＬＳとは異なる観点で）導入された正則化に基づく各種の方法、特に、古くからあるＬ２正則化（リッジ回帰）に加え、多数の説明変数候補の中から必要な説明変数を自動選択することのできるＬ１正則化手法であるLassoなどの方法、Lassoとは異なる視点でベイズ推論的な手法で説明変数の自動選択を行うRVM（適合ベクトルマシン）に基づく回帰手法（予測手法）であるＲＶＲ（適合ベクトル回帰）などの方法、最適化を用いて識別問題（分類問題）をロバスト化したＳＶＭ（サポートベクトルマシン）を回帰問題（予測問題）に転用したＳＶＲ（サポートベクトル回帰）などの方法、従来の重回帰などのような線形回帰を非線形に拡張し、さらに確率分布を考えてノンパラメトリック回帰問題に帰着させ、点予測ではなく予測分布を出力するガウス過程回帰などの方法、あるいは、様々な予測サブシステムを組み合わせて予測性能を向上させるバギングやブースティングと呼ばれる各種の方法（ランダムフォレスト、アダブースト、ＸＧブーストなど）など、多様なアドバンストな時系列データの予測手法も用いられる様になってきている。
　このように、アドバンストなＡＩに基づく予測手法が次々と開発され、従来から用いられる方法と比較して各段に予測精度が向上するという報告も数多くなされる様になり、ＡＩ関連技術に関する期待が高まっている。
　一方で、ＡＩ関連技術が実際の現実的なシステム（インフラシステムやプロセス系のプラント）において、必ずしも容易に展開できない理由の一つとして、統計的およびＡＩ的な解析手法が持つ本質的な性質であるブラックボックス的なアプローチの欠点として指摘されている、「結果に対する説明性の欠如」の問題がある。そして、この「結果に対する説明性の欠如」の問題に対する意識が急速に高まっており、この問題を解決すべく、最近では「説明可能ＡＩ（ＸＡＩ：eXplainable AI)の概念が米国ＤＡＲＰＡにより提唱され、注目が集まっている。
日本国特許第６２６１９６０号公報
　インフラ系のシステムでは、物理的、化学的な知見に立脚したホワイトボックス的なモデルを用いる事も多いが、これは、「Accountability（説明責任）」が求められることに関連していると考えられる。実際、ホワイトボックス的なモデルによる予測精度がブラックボックス的なモデルによる予測精度よりも必ずしも高いわけではなく、また、ホワイトボックス的なモデルが必ずしも、実際の物理現象を正確かつ忠実に再現しているとも限らない。それにも関わらず、ホワイトボックスモデルが時として好まれ、利用される理由は、ホワイトボックスモデルに含まれる「パラメータ」には物理的および化学的な意味があるため「合理的な説明がしやすく、他者に納得してもらいやすい」ことに加え、このような物理的および化学的な意味を持つことにより、精度が十分でない場合に、どのように調整すればよいかという方針や指針をたてやすい、すなわち、調整しやすい、という利点があるためであると考えられる。実際、物理的および化学的に意味を持つパラメータは、パラメータの値自体に物理的および化学的意味があるため、理解がしやすく、そのパラメータの値の取りうる範囲を想定することや、値の増減が結果に与える影響を物理的および化学的知見に基づいて考えることができる。
　上記のように、「合理的な説明のしやすさと調整のしやすさ」がホワイトボックスモデルを用いることの大きな動機付けになっている可能性がある。もしＡＩ的な手法によるブラックボックスモデルにおいても「合理的な説明のしやすさと調整のしやすさ」という要素を付加することができれば、その応用範囲は大きく広がると考えられる。
　本発明の実施形態は、上記事情を鑑みて成されたものであって、予測結果を合理的に説明するとともに、容易に予測結果を調整可能とするモジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラムを提供することを目的とする。
　実施形態によるモジュラー型時系列データ予測装置は、複数のプロセス変数を所定の周期で計測する複数のプロセスセンサを有するシステム又はプロセスに適用される装置であって、複数の前記プロセス変数の時系列データを所定の周期で収集し保存するとともに、複数の前記プロセス変数の中から予測対象となる少なくとも一つの出力変数を選択する出力変数データ選択部と、複数の前記プロセス変数の中から複数の入力変数の候補を選択する入力変数データ選択部と、を含むデータ収集保存部と、前記時系列データから抽出された前記出力変数と複数の前記入力変数との同定用データを用いて、１入力１出力のペア毎にペアワイズ予測モデルのパラメータを同定して複数の前記ペアワイズ予測モデルを定義するペアワイズ予測モデル同定部と、複数の前記ペアワイズ予測モデルから出力されるペアワイズ予測値の合成法を定義する予測モデル合成法定義部と、時間の進行方向に所定の周期又はリアルタイムで前記時系列データから抽出された複数の前記入力変数の予測用データを、複数の前記ペアワイズ予測モデルに入力して、複数の前記入力変数のそれぞれに対応する前記ペアワイズ予測値を演算するペアワイズ出力変数予測部と、前記合成法により複数の前記ペアワイズ予測値を合成して前記出力変数の予測値を演算する合成出力変数予測部と、を有する。
図１は、一実施形態のモジュラー型時系列データ予測装置を適用した雨量流入予測システムを概略的に示す図である。 図２は、一実施形態のモジュラー型時系列データ予測装置の一構成例を概略的に示す図である。 図３は、幹線水位データと流入量データとの関係の一例を概略的に示す散布図である。 図４は、幹線水位データと流入量データとの関係の一例を概略的に示す散布図である。 図５は、一実施形態のモジュラー型の時系列予測装置におけるペアワイズ予測モデルの第１実施例について説明するための図である。 図６は、一実施形態のモジュラー型の時系列予測装置におけるペアワイズ予測モデルの第２実施例について説明するための図である。 図７は、第２実施例のペアワイズ予測モデルの効果を説明するための図である。 図８は、第２実施例のペアワイズ予測モデルの効果を説明するための図である。 図９は、一実施形態のモジュラー型の時系列予測装置におけるペアワイズ予測モデルの第３実施例について説明するための図である。 図１０は、第３実施例のペアワイズ予測モデルの効果の一例を説明するための図である。 図１１は、第３実施例のペアワイズ予測モデルの効果の一例を説明するための図である。
実施形態
　本実施形態では、調整可能な時系列デー予測を行うための基本的手段として、モジュラー型の時系列データ予測装置について説明する。本実施形態のモジュラー型の時系列予測装置は、予測モデル（システム）を部分モデル（サブシステム）に分解可能にして、各部分モデルを容易に着脱可能にするものである。
　モジュラー型の時系列データ予測装置によれば、例えば、複数の監視変数中のいずれか一つの変数のデータが欠測していたり、複数の監視変数にアウトライア（異常データ）が含まれていたり、監視変数がノイズ等を多量に含む不安定な計測データである場合、その影響が全体に波及することを回避できる。また、モジュラー型の時系列データ予測装置によれば、例えば、説明変数が必要であるか否（不要）かを適宜判断して、部分モデルを容易に脱着したり、予測モデルを部分的に調整したりすることができる。
　このような利点は、モジュラー型の構成によりモデルが「分解可能」であるという事によって得られる利点である。予測モデルにおいてモジュラー型の構成をとる利点は、単に分解可能ということによって得られる利点だけでなく、これに加えて、このモジュラー型の構成をとることで、物理法則に矛盾せずに合理的で納得できるモデルを構築することが可能になり、さらに、このことにより調整を容易にすることができる点である。
　予測モデルにおいてモジュラー型の構成をとる利点について以下に説明する。
　一般に、ブラックボックス的な手法が説明可能（Explainable）でない理由は、入力および出力という外部情報だけしかわからず、内部構造が見えない不可視である事が原因であると考えられている。すなわち、ブラックボックス的な手法は、単に内部構造が見えないだけでなく、そもそも、外部の入出力関係からだけでは、内部構造を一意に決めることができず、同じ入出力関係を与える複数（無数）の内部構造がありうる。このため、内部構造の解釈が難しくなり、場合によっては物理法則に矛盾する不合理なモデルになってしまうことがある。これは、モデルの一意性の問題として、制御理論分野の中に含まれるであるシステム同定分野等で知られた事項である。このようなモデルの一意性の概念は、例えば、「ある入力データと出力データのペアに対して、その入出力関係を表すモデルを唯一に決定することができるか」という問題であり、モデルの可同定問題と呼ばれている。
　可同定問題は、モデルの構造（モデルを表現する数式の形）を予め限定した場合、例えば、モデルを線形伝達関数モデルや時系列モデルに限定した場合、そのモデルに含まれるパラメータ（例えば伝達関数に含まれる係数）を一意に決めることができるかという問題となる。この場合、可同定問題の中でも特にパラメータ可同定性の問題と呼ばれる。このパラメータ可同定性の問題は、現実のデータを扱う問題でも頻繁に現れる問題であり、以下の３つのケースは特に遭遇することが多い問題である。
　一つ目は、先に述べた多重共線性の問題であり、これは、入力変数（説明変数）間に強い相関を持つ場合の問題である。これは、例えば、雨水ポンプ場への雨水流入量を予測する場合に、入力変数としてレーダ雨量情報を用いて、近接する異なる２つのメッシュの降雨量を入力変数とするような場合などに典型的に現れる問題であり、現実のデータでは、説明変数間に強い相関を持つことは極めて多い。
　二つ目は、制御理論の分野でよく知られている、持続的励振条件（Persistent Exciting Condition、PE条件）と呼ばれる条件が成立しない場合の問題であり、これは、モデル化の対象プロセスに対して、十分な周波数成分を含まない入力データが適用されている場合に生じる問題である。極端ではあるが、現実にも認められる例として、入力データが一定値である場合、通常出力データも一定値になるが、このような一定値の入出力関係を表す伝達関数モデルは無数に存在するため、係数を一意に決めることができない。例えば、y(t)=b1×u(t-1)+b2×u(t-2)と表される簡単な伝達関数モデル（時系列モデル）を考える場合、u(t-1)=u(t-2)=ｋ（一定値）であるとb1とb2を一意に決めることはできない。この概念を一般化したものが持続的励振条件というパラメータ可同定条件であるが、このような条件が成立しない場合も多い。
　三つ目は、モデルの構造の複雑さと比較して、データの得られる量が少なくデータの質が良くない場合である。典型的な例として、ある対象が非線形的な挙動を示すことが物理的な観点から予めわかっているが、得られるデータが少なかったり、ある特定の運用条件の付近のデータしか得られなかったりすることにより、非線形の係数を同定できない様な場合である。例えば、y(t)=c1×u(t)+c2×u(t)^2と書かれることがわかっているが、データとして得られるu(t)の値が小さい値の場合、u(t)^2<<u(t)となるため、c2を同定することが困難になる様な場合である。
　上記三つの例は、パラメータの可同定性が欠落する典型的なケースであるが、現実の問題では、このような状況に遭遇することが極めて多い。モデル化の対象が、比較的小さな機械系システムおよび電気系システムの場合は、可同定性を向上させるための実験を行い、可同定性を確保した上で、モデルを構築することが常道である。しかし、上下水道、雨水排水、電力システム、などの巨大な社会インフラシステムの場合、モデル化の対象となるデータが自然現象から得られるデータ（降雨、日射量、風速など）であったり、対象プラントが稼働している状態のデータであったりすることが多く、パラメータの可同定性を向上させるための実験などを行うことができず、モデル化に利用できるデータは与えられたものであることが多いため、可同定性の欠落を回避することが本質的に難しい場合が多い。この結果、社会インフラシステムの場合には、可同定性が欠落した状態でモデルを構築せざるを得ない。
　以下、可同定性の欠落による問題の一例として「多重共線性」について、従来の技術と本実施形態のモジュラー型時系列データ予測装置とを対比する。
　「多重共線性」とは、複数の説明変数（入力変数）間に相関関係（一次従属関係）があることにより、回帰モデルのパラメータを一意に決めることができなくなるという問題である。これを最も簡単な例を用いて説明する。
　一つの出力yが二つの入力ｕ１、ｕ２の重回帰で表される以下の重回帰モデルを考える。
　　y=a1×u1＋a2×u2　　（１）
　(1)式は最も簡単な重回帰モデルであり、ａ１とａ２とは回帰係数と呼ばれるパラメータである。ここで、物理的な具体的イメージを持つために、例えば、出力ｙを雨水ポンプ場への雨水流入量、入力ｕ１および入力ｕ２を各々異なる箇所の降雨量であることを想定して考えてみる。この場合、降雨量が増加すれば流入量が増加することは、物理的に考えて自明であるから、ａ１とａ２とは正の値を持たなければならないことは容易にわかる。ここで、降雨量ｕ１と降雨量ｕ２とに従属関係のある多重共線性の問題がある場合を想定する。
　この多重共線性問題は、実際の現象においても頻繁に見られることであり、この場合、例えば、入力ｕ１と入力ｕ２とをレーダ雨量データの隣接メッシュの降雨量であることを想定してみれば、この問題が特殊な問題ではなく、むしろ自然な問題であることが容易にわかる。例えば、国土交通省が配信するＸＲＡＩＮと呼ばれる気象レーダでは、１メッシュのサイズが２５０ｍ×２５０ｍであるため、隣接メッシュの降雨量とは互いに数百メートル離れた箇所の降雨量を意味しており、その降雨量がほぼ同じであることは物理的に考えて極めて自然な事である。２つの降雨量データがレーダ雨量の隣接メッシュの降雨量ではなくて、２か所の地上雨量データであったとしても、通常雨水排水区の大きさが数キロ四方程度以内であることを考えると比較的近い値を持つことは容易に想像できる。従って、このような場合ｕ１≒ｕ２となっていることは、自然な事であり、これは多重共線性問題そのものである。
　そこで、ここでは簡単のため、ｕ１＝ｕ２という完全な多重共線性の状態にあった場合を想定する。そして、例えば、２０ｍｍ／ｈの強度の雨量に対して、１０ｍ３／ｓの雨水流入があるとする。今、ｕ１＝ｕ２であるので、ｕ：＝ｕ１＝ｕ２と定義すると、（１）式は、以下の（２）式の様に書き換えられる。
　　y=(a1+a2)×u　　（２）
　上記(2)式は単回帰の形であり、ｕ：＝ｕ１＝ｕ２であるから、ａ１＋ａ２＝０．５となる。しかし、ａ１＋ａ２＝０．５という制約条件を満足してさえいれば、ａ１とａ２とがいかなる値であっても、（２）式は全く同じ出力結果を出すため、ａ１とａ２との各値はａ１＋ａ２＝０．５を満たす範囲で任意の値を持つことができる。そのため、例えば、ａ１＝３．５、ａ２＝－３でも、ａ１＝１００、ａ２＝－９９．５となっても問題がなく、全く同じ入出力関係を与えるため、各々のパラメータ値の相違を入出力データから識別することができなくなる。
　一方、ａ１とａ２とは、各雨量に対する係数であるから、本来、すくなくとも正の値を持つべきであり、２地点の降雨量が同じであるなら、その値はほぼ等しくなるべきであることも物理的に考えれば妥当である。しかし、このように、入力変数間の相関関係（一次従属関係）により、パラメータ値を一意に決めることができず、本来のパラメータ値を正しく同定できなくなる事が多重共線性の問題の本質である。
　このような多重共線性などのパラメータ可同定性の欠如により、入出力関係から一意にパラメータを決めることができなくなる問題は、システム同定分野や統計分野（統計的機械学習分野≒AI分野）では、数学的な「悪条件問題」として扱われ、この悪条件を回避するための各種のアルゴリズムが開発されている。先に述べた正則化やＰＬＳなどの方法は、このような悪条件回避のための代表的な技術的手段（テクニック）であり、アドバンストなＡＩ手法の多くには、正則化や（ＰＬＳで用いられる）次元圧縮などの手法(テクニック)が採用されている。
　しかし、例えば、正則化やＰＬＳなどのテクニックを駆使したアドバンストなアルゴリズムを（１）式の問題に適用したとしても、上記の問題に関して、ａ１＞０、ａ２＞０、かつａ１≒ａ２となるような値を推定するとは限らず、むしろ、そのような値にならない可能性の方が高い。これは、正則化や次元圧縮などは、多重共線性（やこれを一般化した可同定性の欠如）の問題を、数学的な悪条件の問題として捉え、これを回避するためのある種の“テクニック”であり、悪条件を回避してパラメータを推定（同定）できるようにしているだけであるからである。すなわち、多重共線性を持つ重回帰問題に対して古典的な最小２乗法を適用すると、「答えを求めることができない＝パラメータ値を求めることができない」という問題になるため、悪条件を回避するというテクニック（正則化や次元圧縮）をアルゴリズムに導入することで、「答えを求めることができる＝パラメータ値を同定できる」問題にしているだけであり、同定した値が物理的に合理的に納得できるものであるか否かを直接考慮しているわけではないからである。
　実際、例えば、正則化のテクニックでは、悪条件を回避するために、パラメータを同定する際に、予測誤差を評価するだけでなく、「パラメータの値自身が大きくなりすぎない」というパラメータに関するノルム（大きさを図る指標）を小さくする、という評価指標を加えることで悪条件を回避しており、これはパラメータの物理的な意味づけや解釈を行うこととは直接関係していない。
　一方、このような問題に対して、合理的に納得のできるパラメータ値を得る方法について以下に説明する。はじめに、（１）式のモデルではなく、（３）式と（４）式との単回帰モデルを考える。
　　y=a1´×u1　　（３）
　　y=a2´×u2　　（４）
　この時、ｕ１＝ｕ２であり、ｕ１＝ｕ２＝２０ｍｍ／ｈに対して、ｙ＝１０ｍ３／ｓを出力するモデルは、単位系をそのままで考えると、ａ１´＝ａ２´＝０．５であることは明らかである。従って、このような関係を持つ場合に、（３）式と（４）式と単回帰を適用すると、ａ１´＝ａ２´＝０．５に近い値がパラメータとして同定（学習）されることになる。この０．５という値は、（回帰モデルなので物理的な意味は明確ではないが）合理的で納得できる値である。
　次に、（３）式と（４）式とを用いて、（１）式の重回帰モデルを構築することを考える。
　今回の例では、入力変数が２か所の降雨量となり、その和をとっているため、入力となる降雨量が倍になるので、各降雨量に対する係数のａ１とａ２とは、一つの降雨量を入力とした単回帰モデルの係数の半分となり、ａ１＝ａ２＝０．２５となると考える事が物理的に妥当な納得できるパラメータ値であると考えられる。
　このようにして、（３）式と（４）式とでｕ１とｙとの関係、および、ｕ２とｙとの関係を単回帰で各々求めて、２つの予測モデルを構築する。そして、次に２つの予測モデルの出力ｙの平均を計算する。このように、構成的（Constructive）な手順を踏むと、ａ１とａ２との値は、各々０．２５と、物理的に妥当と考えられる納得のできる値として推定され、結果的に、（１）式の形の重回帰式の係数を定めることができる。これは、部分モデル（サブシステム）を結合していくモジュラー型のアプローチそのものである。
　上記重回帰モデルにより係数を定めることができる理由は、一般に、入出力関係を表すデータとして同じデータを用いる場合、モデル構造が簡単（≒モデルに含まれるパラメータ数が少ない）なモデルの方が、可同定性が高くなる事（＝パラメータ値を一意に同定可能）が知られているためである。本実施形態のモジュラー型時系列データ予測装置は、上記の特性を利用している。すなわち、同じ入出力データであったとしても、その可同定性はモデルの構造やパラメータ数に依存し、パラメータ数が多いモデルでは可同定でなくなる場合であっても、パラメータ数が少なければ可同定になる場合がある、ということを利用している。
　上記の例でいえば、本実施形態のモジュラー型時系列データ予測装置では、２入力の重回帰モデルでは多重共線性により可同定でなくなるが、これを１入力の単回帰の組み合わせと考えることにより、多重共線性の問題を本質的に回避して、可同定なモデルを積み上げることによって合理的なモデル構築を可能にしている。
　次に、本実施形態のモジュラー型時系列データ予測装置と、いわゆる「Explainable　AI（ＸＡＩ）」との関係および相違について説明する。
　いわゆるＸＡＩでは、全体を一つの大きなシステムとして捉えてブラックボックス的にモデルを構築して、後から、そのモデルを適切に解釈可能な形に分解して説明を加えようとするトップダウン型のアプローチ（あるいは帰納的なアプローチ）である。これに対し、本実施形態のモジュラー型時系列データ予測装置の考え方は、個別には比較的簡単に理解できる合理的な説明が可能なモデル（サブシステム）を統合していく事により、全体の予測モデルを構築していくボトムアップ型のアプローチ（あるいは演繹的なアプローチ）である。
　上記のようなボトムアップ型のアプローチでは、サブシステム間の相互の干渉を予め考慮することができないため、予測精度面に注目した場合には、トップダウン的（ＸＡＩ的）なアプローチに劣る可能性を否定できない。しかしながら、ボトムアップ型のアプローチは、合理的な説明性という意味では、先に説明したようにパラメータの可同定性を維持したモデル構築が可能になり、トップダウン的なアプローチよりも優れる可能性が高い。
　また、ボトムアップ型のアプローチは、ホワイトボックス的なモデリングのアプローチとも類似している。ホワイトボックスモデリングは、基本的に演繹的、構成的なモデリング手法である。ホワイトボックスモデリングでは、物理的に意味のあるパラメータを用いて構成的にモデルを構築するのに対し、本実施形態のモジュラー型アプローチでは、物理的な意味を持たせることはできないものの、個々のパラメータを入出力と直接関係づけることにより解釈が容易なパラメータを同定し、それを用いて構成的にモデル構築している。
　本実施形態のモジュラー型アプローチがＸＡＩ的なアプローチよりも合理的な説明性において優れていることは以下の事からも明らかである。例えば、先の例において、もしトップダウン的（ＸＡＩ的）なアプローチを行おうとすると、正則化やＰＬＳなどの方法を使って係数パラメータのａ１やａ２を求め、それに対して何等かの解釈を加える処理を施して結果を「Explainable（説明可能）」にする方法の確立を目ざす。しかしながら、上記の簡単な例からわかる様に、ａ１とａ２が０．２５という値からかけ離れた値で同定されていた場合（実際、正則化やＰＬＳを用いると、そのようになる可能性がある）には、いかなる方法を用いて「説明」を試みても、合理的な説明が原理的に困難であることは容易にわかる。例えば、ＰＬＳや正則化を用いてパラメータ値を求めることができても、ａ１かａ２のいずれかが負の値（例：ａ１＝３．５、ａ２＝－３）となっている場合には、合理的に納得のできる説明は不可能である。
　つまり、本来合理的に納得できる値である0.25という値を求められていなかったとしても、可同定性が欠如している場合には、ある制約さえ満たしていればどのような内部パラメータ値であろうが予測精度だけに着目する限りは差異が見られないため、「精度の良いモデル」と見なされる。しかしながら、精度が良いモデルが、必ずしもそのモデルのパラメータ値の妥当性を補償する合理的で妥当なモデルを意味しないため、トップダウン的（ＸＡＩ的）なアプローチでは原理的に合理的な説明が困難になってしまう。
　なお、ＸＡＩ的なアプローチでは、ａ１とａ２とのそれぞれに対する個別の意味付けは困難であるが、ａ１＋ａ２に対する意味づけは可能であるため、説明不可能であるというわけではない。ただし、ＸＡＩ的なアプローチでは、内部に含まれる全てのパラメータに対して意味付けを行うことは原理的に困難になる可能性が高く、何等かの「特徴量」を抽出して、それに対して意味付けを求めることになる。
　従って、「合理的な説明性の高いモデル」を構築することを主目的にするのであれば、そもそも、正則化やＰＬＳなどの数学的に悪条件を回避する様なアルゴリズムに頼る前に、可同定性が欠如する様な状態を回避してモデルの内部構造（内部パラメータ）の一意性を保証する様にモデルを構築しておき、合理的な説明性を担保しながらモデルを統合していくボトムアップ的、構成的、演繹的なモジュラー型アプローチの方が優れていると可能性が高い。
　また、このように合理的説明性を担保したモデルを構築できることに加え、最初に述べた分離可能な構造とすることで、モデルの調整も容易にすることができ、これらを合わせて合理的説明可能かつ調整可能なモデル構築を行うことが可能になる。
　本実施形態のモジュラー型時系列データ予測装置は、以上の洞察に基づいて成されたものである。
　以下、実施形態のモジュラー型時系列データ予測装置について、図面を参照して詳細に説明する。
　図１は、一実施形態のモジュラー型時系列データ予測装置を適用した雨量流入予測システムを概略的に示す図である。
　以下では、本実施形態のモジュラー型時系列データ予測装置を搭載したシステムの実施イメージとその効果をより明確にするために、雨水流入予測システムを対象として説明する。なお、本実施形態のモジュラー型時系列データ予測装置が適用される対象プロセスは雨水流入予測システムに限定されるものではなく、２種類以上の計測項目の時系列データが存在する任意のプロセスに対して適用することができる。モジュラー型時系列データ予測装置は、例えば、鉄鋼プラント、石油化学プラント、食品プラント、製薬プロセス、半導体製造プロセス、発電プラント、交通監視設備、空調監視設備、などのプロセスに適用することができる。
　雨量流入予測システムの対象プロセス１は、流量計１１と、幹線流量計１２と、幹線水位計１３１-１３Ｋと、地上雨量計１４１-１４Ｍと、レーダ雨量計１５と、ポンプ井水位計１６１、１６２と、流入渠１７と、雨水ポンプ井１８と、雨水ポンプ１９と、流入ゲート１１０と、を含む都市雨水排水プロセス１である。レーダ雨量計１５は、Ｑ×Ｐメッシュのメッシュ毎のレーダ雨量計１５１１～１５ＱＰを含む。
　都市雨水排水プロセス１の各種センサ１１－１５、１６１、１６２は、所定の周期（例えば、３０秒、６０秒）で、プロセスの状態を表す量や運転操作に関わる量などの計測を行う。各種センサ１１－１５、１６１、１６２で計測された値は、モジュラー型時系列データ予測装置で収集および保存される。
　流入渠１７は、下水管やポンプ場から送られてきた雨水が流れ込む渠である。
　雨水ポンプ井１８は、流入渠１７から流入した雨水が貯められる貯水槽である。雨水ポンプ井１８に流入する雨水は、沈砂池で一緒に流れてきた砂などが予め取り除かれていてもよい。
　雨水ポンプ１９は、雨水ポンプ井１８に貯められた雨水を強制的に河川等へ送り出すポンプである。
　流入ゲート１１０は、流入渠１７と雨水ポンプ井１８との間の流路に設けられ、開閉することにより流入渠１７から雨水ポンプ井１８に流入する雨水の量を調整するように動作する。流入ゲート１１０の前後には水位計１６１、１６２が設置され、例えば、水位計１６１、１６２により測定された値に応じて流入ゲート１１０の動作が制御される。
　図２は、一実施形態のモジュラー型時系列データ予測装置の一構成例を概略的に示す図である。
　本実施形態のモジュラー型時系列データ予測装置は、データ収集保存部２と、データ抽出部３と、ペアワイズ予測モデル同定部４と、予測モデル合成法定義部５と、ペアワイズ出力変数予測部６と、合成出力変数予測部７と、予測誤差評価部８と、ペアワイズ予測モデル修正部９と、出力予測結果観測部１０と、を備えている。
　本実施形態のモジュラー型時系列データ予測装置の構成は、例えば、ＣＰＵやＭＰＵなどのプロセッサを少なくとも１つと、プロセッサにより実行されるプログラムが格納されるメモリとを備えた演算装置であって、種々の機能をソフトウェアにより、若しくは、ソフトウェアとハードウエアとの組み合わせにより実現することが出来る。
　データ収集保存部２は、各種センサ１１－１５、１６１、１６２で計測されたプロセス変数の値を収集し保存する。データ収集保存部２は、出力変数データ選択部２１と、入力変数データ選択部２２とを備えている。
　出力変数データ選択部２１は、各種センサ１１－１５、１６１、１６２で計測された値の少なくともいずれかの量を出力変数として選択する。
　入力変数データ選択部２２は、各種センサ１１－１５、１６１、１６２で計測された値から、複数の入力変数として選択する。なお複数の入力変数の数は、複数のプロセス変数の数以下であればよい。
　データ抽出部３は、データ収集保存部２に保存された時系列データの中から所定のデータを抽出する。データ抽出部３は、オフライン予測モデル同定用データ抽出部３１と、オンライン予測用データ抽出部３２と、評価用データ抽出部３３とを備えている。
　オフライン予測モデル同定用データ抽出部３１は、データ収集保存部２から、所定の同定周期若しくは外部からの要求により、予測モデル同定用のデータとして各種プロセスセンサ１１－１５、１６１、１６２で計測された（当該時刻よりも過去の）所定の期間のデータを抽出する。
　オンライン予測用データ抽出部３２は、データ収集保存部２から、所定の予測周期で、各種プロセスセンサ１１－１５、１６１、１６２の計測値からオンラインの予測に必要となる予測用データを、リアルタイムで抽出する。
　評価用データ抽出部３３は、データ収集保存部２から、所定の評価周期もしくは外部からの要求により、予測誤差評価に必要となる各種プロセスセンサ１１－１５、１６１、１６２の計測値を、リアルタイムで抽出する。
　ペアワイズ予測モデル同定部４は、オフライン予測モデル同定用データ抽出部３１で抽出された計測値のデータを用いて、出力変数データ選択部２１で選択された出力変数および入力変数データ選択部２２で選択された入力変数を変数毎にペアワイズとし、後述する方法により予測モデルを同定する。
　予測モデル合成法定義部５は、ペアワイズ予測モデル同定部４で同定したペアワイズの予測モデルの出力値（ペアワイズ予測値）から、出力変数（出力変数が複数の場合は出力変数毎）を合成する方法を定義する。
　ペアワイズ出力変数予測部６は、ペアワイズ予測モデル同定部４で同定したペアワイズの予測モデルに対して、時間の進行方向における所定の監視周期で、オンライン予測用データ抽出部３２から抽出した入力変数各々に対応する予測用データを入力して、入力変数毎に出力変数（出力変数が複数の場合は出力変数毎）のペアワイズの予測を行う。
　合成出力変数予測部７は、ペアワイズ出力変数予測部６の出力であるペアワイズの予測値に対して、予測モデル合成法定義部５で定義した方法に従って合成した出力変数の予測値を出力する。
　予測誤差評価部８は、ペアワイズ出力変数予測部６の出力であるペアワイズの予測値と、合成出力変数予測部７の出力であるペアワイズの予測値とを合成した予測値を保存し、所定の周期もしくは所定のタイミングで、評価用データ抽出部３３により抽出した出力変数の実績値と保存したペアワイズおよび合成予測値との誤差を、各々評価する。
　ペアワイズ予測モデル修正部９は、予測誤差評価部８によって評価した予測誤差に基づいて、ペアワイズ予測モデルの調整の必要性の有無を判断し、必要なペアワイズ予測モデルの修正を行う。ペアワイズ予測モデル修正部９は、例えば、ペアワイズ予測モデルのパラメータ調整、もしくは、ペアワイズ予測モデルの削除、もしくは、ペアワイズ予測モデルの統合、を行うことができる。ペアワイズ予測モデル修正部９は、修正したペアワイズ予測モデルをペアワイズ出力変数予測部６に供給する。
　出力予測結果観測部１０は、例えばモニタなどの表示手段と、ユーザインタフェースとを備え、合成出力変数予測部７による合成出力変数の予測値を少なくとも観測者に提示することができる。出力予測結果観測部１０は、例えば、パーソナルコンピュータ、タブレット端末や、スマートフォンなどの携帯端末であってもよい。出力予測結果観測部１０は、必要に応じて、ペアワイズ出力変数予測部６によるペアワイズ予測モデルの複数の出力と、それから合成される予測分布情報の少なくとも一方を、観測者に更に提示してもよい。
　次に、本実施形態のモジュラー型時系列データ予測装置の動作について説明する。
　例えば都市雨水排水プロセス１では、各種プロセスセンサ１１－１５、１６１、１６２によって、所定の周期でプロセスの情報が計測され、データ収集保存部２に計測値（プロセス変数の値）が供給される。データ収集保存部２は、収集した計測値を、予め決められたフォーマットに従って、時系列データとして保存する。
　出力変数データ選択部２１は、各種プロセスセンサ１１－１５、１６１、１６２で計測されたプロセス変数の中から予測の対象となる変数を選択して、選択した変数を出力変数と設定する。
　例えば都市雨水排水プロセス１では、通常雨水ポンプ井１８への流入量を予測の対象とする場合には、出力変数データ選択部２１は、流量計（雨水ポンプ井流入量計）１１で計測されている雨水ポンプ井流入量が予測対象の出力変数として選択する。
　また、流量計（雨水ポンプ井流入量計）１１が設置されていない様な場合には、雨水ポンプ井水位計１６１、１６２の計測値と、雨水ポンプ１９および流入ゲート１１０の運用状態とから、雨水ポンプ井１８への流入量を計算によって求めて、これを予測対象として選択しても良い。この場合、出力変数データ選択部２１は、流入量の計算値を計算して保存しておき、流入量の計算値を出力変数として選択する。
　また、例えば都市雨水排水プロセス１の複数の項目を予測する場合には、出力変数データ選択部２１は、複数の項目を出力変数として選択してもよい。都市雨水排水プロセス１においても、例えば、雨水ポンプ井１８に複数の箇所から雨水の流入があるような場合は、出力変数データ選択部２１は、複数の流入箇所各々における流入量を出力変数として選択してもよい。
　入力変数データ選択部２２は、出力変数データ選択部２１で選択した出力変数に影響を与える可能性のある測定値を入力変数の候補として選択する。影響を与えるか否かについて事前情報が全く無い場合には、入力変数データ選択部２２は、全てのプロセス変数を機械的に入力変数としておいても良い。入力変数データ選択部２２は、出力変数として選択した変数を入力変数としても選択して構わない。例えば、時系列のデータ解析では、出力変数として選択した変数自身の過去の値から未来の値を予測する所謂自己回帰を検討することも多くある。このような場合には、流量計（雨水ポンプ井流入量計）１１による雨水ポンプ井流入量は、出力変数と入力変数との両方として選択されても構わない。
　また、出力変数データ選択部２１で選択した出力変数が複数ある場合、入力変数データ選択部２２は、ある出力変数に対して別の出力変数を入力変数として選択してもよい。このようにすることで、複数の出力変数同士の相関を考慮することが可能になる。一方、複数の出力変数が異なる流入箇所の流入量であるときには、プロセスの物理的な構造上、出力変数同士の関連が無い場合が多い。この場合、入力変数データ選択部２２は、ある出力変数に対する入力変数として別の出力変数を除外しておくこともできる。
　入力変数データ選択部２２は、その他の変数に対しても、例えば、物理的な管渠の接続関係などの情報から明らかに出力変数に無関係な変数や因果関係（入力（原因）と出力（結果）との関係）が物理的に成立しない変数を予め入力変数から除外しておくこともできる。
　例えば、雨水ポンプ井水位計１６１、１６２により計測される雨水ポンプ井水位は、出力変数として選択した流量計（雨水ポンプ井流入量計）１１による雨水ポンプ井流入量の影響を受けて変化することは物理的に考えて明らかであるから、雨水ポンプ井水位を入力変数から除外しておく。
　入力変数データ選択部２２は、物理的な構造の関係が明確でない変数について予め出力変数データ選択部２１で選択した出力変数との相関解析を予備的に実施し、相関係数の絶対値が所定の値（例えば０．５）未満の変数を、予め入力変数から除外しておいてもよい。さらに、入力変数データ選択部２２は、相関解析を実施する際に、解析対象の変数が測定された時間をずらしながら相関解析を行い、最も相関の高くなる時のずらした時間Ｌを求め、時間Ｌの符号情報から、因果関係を満たさない符号を持つ変数を入力変数から予め除外しておいてもよい。なお、先の雨水ポンプ井水位の変数は、このような解析によっても入力変数から除外される可能性が高い。
　本実施形態のモジュラー型時系列データ予測装置では、入力変数データ選択部２２は、雨水水位計１６１、１６２により計測された雨水ポンプ井水位のみを除いた流量計（雨水ポンプ井流入量計）１１による雨水ポンプ井流入量と、幹線流量計１２により計測された幹線流入量と、幹線水位計１３１-１３Ｋにより計測された幹線水位と、地上雨量計１４１-１４Ｍにより計測された地上雨量と、Ｑ×Ｐメッシュのメッシュ１５１１-１５ＱＰ各々におけるレーダ雨量とを入力変数として選択するものとする。
　なお、本実施形態のモジュラー型時系列データ予測装置において、例えば、Ｋ個の幹線水位計１３１～１３Ｋは、最近、下水管渠内の水位情報を正確に把握しようとする行政の動向などから、Ｋの値が１０～１００程度の比較的多数の幹線水位計が設置されることも多く、また、気象レーダは国土交通省が全国にＸＲＡＩＮという気象レーダを設置していることから、多くの場所で取り入れることができ、対象とする雨水排水区の大きさにもよるが、Ｑ×Ｐメッシュのメッシュサイズは１００程度になることも珍しくないため、実際に計測される時系列データの項目数（変数の数）は、かなり多くなる。また、本実施形態のモジュラー型時系列データ予測装置が適用される都市雨水排水プロセス１に限らず、計測変数の項目数が数百～数千に及ぶことは一般的である。
　本実施形態のモジュラー型時系列データ予測装置では、入力変数は、例えば、雨水ポンプ井水位計１６１、１６２による計測値のみを除いた流量計（雨水ポンプ井流入量計）１１による雨水ポンプ井流入量の計測値と、幹線流量計１２による幹線流入量の計測値と、幹線水位計１３１～１３ＫによるＫ個の幹線水位の計測値と、地上雨量計１４１～１４ＭによるＭ個の地上雨量の計測値と、Ｑ×Ｐメッシュ１５１１～１５ＱＰの各メッシュのレーダ雨量の計測値と、である。以下では、入力変数の数をｐ（正の整数）として説明する。
　オフライン予測モデル同定用データ抽出部３１は、予測モデルを構築する人が指定した所定期間のデータを用いて要求されたタイミングで、先に出力変数データ選択部２１と入力変数データ選択部２２との該当する所定期間分の時系列データを取得する。オフライン予測モデル同定用データ抽出部３１は、予測モデルの構築を定期的に行い、モデルを定期的に更新したい場合には、時間の進行方向において所定の周期ＴＬで所定期間の入出力変数の時系列データを定期的に抽出しても良い。このようにして抽出されたデータセットをｕｋ、ｋ＝１、２、…ｐと記載する。ｕｋは列ベクトルであるとし、各行が時系列データの各時刻のサンプルに対応する。同様に出力変数である雨水流入量の時系列データセットｙとしておく、オフライン予測モデル同定用データ抽出部３１は、データセットｕｋとデータセットｙとを所定の周期で抽出する作用を持つ。
　次に、ペアワイズ予測モデル同定部４では、ｕｋとｙとを用いて、（ｕ１、ｙ）、（ｕ２、ｙ）、…（ｕｐ、ｙ）のｐ個のペア毎に予測モデルを同定する。本実施形態のモジュラー型時系列データ予測装置では、出力変数ｙは、流量計（雨水ポンプ井流入量計）１１による雨水ポンプ井流入量のみの一つであるが、複数（ｈ個とする）の出力ｙがある場合には、ｈ個の各々のｙに対して同様の処理を行えば良く、以降の処理は出力変数の数がｌ個に増えた場合でも、ｈ回繰り返して同じ処理を行えばよいだけであるので、一般性を失う事なく出力変数は１個であるとしておく。
　最も典型的な予測モデルは、次式で表せる線形回帰の形をしたモデルを仮定して同定を行うものである。　
　　y(t)=a1×uk(t-L)+a2×uk(t-L-1)+…+an×uk(t-L-n+1)+c　　（５）
　ここで、ａｋ、ｋ＝１、２、…、ｎとｃは同定すべきパラメータである。ｃは入力ｕｋと出力ｙの平均の差を表すバイアスパラメータであるが、予めｕｋとｙから平均値を除去して平均を０としておくことにより、常に０とすることができるので、以下では簡単のためにｃは０としておく。
　上記（５）式において、ｕｋがｙ自身ではない場合を、ディジタル信号処理やシステム同定の分野では、有限インパルス応答モデル（ＦＩＲ）モデルと呼ぶ（入力変数に対する重み付き移動平均モデルと呼ぶこともできる）。一方、ｕｋがｙ自身と一致する場合は自己回帰モデル（ＡＲモデル）と呼ばれる。ここで重要な事は（５）式の線形回帰の形であっても、複数の入力変数による有限インパル応答を足し合わせ、さらに自己回帰を足し合わせた多入力１出力（複数の出力を同時に考える場合は多入力多出力）の線形回帰の形（ＡＲＸ（Autoregressive eXogenous Input）などの形）ではなくて、あくまでも１入力１出力の形の線形回帰モデルになっている点であり、これが「ペアワイズ」の意味である。
　ペアワイズで同定することの第一の利点は、「複数の入力」という概念自身が無くなるため、先に述べた（複数の入力変数間の関係から生じる）多重共線性の問題が決して生じる事はなく、この問題を本質的に回避できることを含め、パラメータの可同定性の観点では多入力１出力のシステムよりも有利であるからである。これを若干補足する。システム同定の分野では、入力信号ｕｋに含まれる周波数成分の量（数）によって同定可能なパラメータの数が変化することが持続的励振条件（ＰＥ（Persistent Exciting）条件）と呼ばれるパラメータの可同定性（＝パラメータを一意に決定可能）条件として知られており、同定（推定）すべきパラメータａｋの数が多くなるほど、入力信号ｕｋは多数の周波数成分を含まなければならない事が知られている。
　一方、本実施形態のモジュラー型時系列データ予測装置において、入力変数の値は、予測モデルの構築者が制御（調整）できるものではなく、所与の観測（計測）情報であるから、どれくらいの周波数成分を含んでいるか否かは、外部条件として与えられるものであり、パラメータの数は少なければ少ないほど可同定性が向上し、逆にパラメータの数が多ければ多いほど可同定性が劣化する可能性が高くなる。より厳密に述べると、パラメータの数が少ない場合に可同定性が悪くなる可能性は０％であり、パラメータ数が少ないものの可同定性の方が高いか、悪くてもパラメータ数が多いものの可同定性と同じである。従って、上記（５）式の１入力１出力のモデルではｎ個のパラメータを含むのに対し、ｐ入力１出力の予測モデルを考えるとｐ×ｎ個のパラメータを含むことになるため、入力数ｐが多くなればなるほどパラメータを一意に同定することが難しくなる可能性が高くなる事がわかる。
　さらに、先に述べた入力変数間の従属関係（相関関係）があると、多重共線性の問題が生じ、同定すべきパラメータ数の多寡に関わらず可同定性が失われ、入力情報は所与のものであるから、いかなる方法を用いても原理的にパラメータ値を一意に決定することが不可能になる。これに対し、１入力１出力の（５）式のモデルでは多重共線性の問題が原理的に生じることはなく、（５）式の形でパラメータが同定できなくなる場合は、入力ｕｋに十分な周波数情報が含まれていない時に限られ、この場合は、次数ｎを調整して同定すべきパラメータ数を調整することで、パラメータａｋを必ず可同定（一意に同定可能）にすることができる。
　例えばｕｋの値が一定値である場合、uk(t-L)=uk(t-L-1)=…=uk(t-L-n+1)となるが、このような場合はｎ＝１として一つのパラメータａ１だけを同定する様にすればパラメータの可同定性を維持することができる。一般的には、所与のｕｋを用いて上記の持続的励振条件を調べることで同定可能なパラメータ数ｎの上限値ｎｍａｘを求めることができるので、ｎ≦ｎｍａｘの範囲内で、例えば、ＡＩＣ（赤池情報量規範）やＢＩＣ（ベイズ情報量規範）、ＭＤＬ基準（最小符号化基準）などの各種の規範および基準を用いたり、交差検証（クロスバリデーション）を行ったりすることで、ｎの適正値を決めて、（５）式の可同定性を常に維持するこができる。
　ペアワイズで同定することの第二の利点は、たとえ、出力変数を除く入力変数が１つの場合であっても、この方法では出力変数を入力とする自己回帰ＡＲモデルと入力変数（出力変数ではない入力変数）を入力とするＦＩＲモデルとの二つのモデルを分離して、同定することによる効果である。このような場合、ペアワイズでない通常の方法では、ＡＲＸモデルなどを用いて、自己回帰成分（ＡＲ成分）と重み付き移動平均成分（ＦＩＲ成分／ＭＡ成分）とを同時に同定することになるが、多くの場合、出力変数と入力変数との関係（相関）よりも出力変数自身の相関（自己相関）の影響の方が強い場合が多く、予測に対する自己回帰成分の影響が極めて強くなる場合が多い。このような場合、入力変数の値が急変した場合にその変化に追従できなくなり、予測が遅れるという現象が生じやすい。
　本実施形態のモジュラー型時系列データ予測装置の場合、例えば、出力変数はポンプ井流入量であり、入力変数として１か所の地上雨量を仮定すると、流入量に対する自己回帰成分の影響の方が地上雨量に対するＦＩＲ成分の影響よりも強くなり、降雨の変化により地上雨量の値が急変した場合でも流入量に変化が現れるまでは予測値の変化が小さく、予測値が実測値に対して遅れるという現象が生じやすい。これに対し、流入量に対する自己回帰モデル（ＡＲモデル）と地上雨量に対する有限インパルス応答モデル（ＦＩＲモデル）とを各々個別に同定するペアワイズの方法では、監視員等のユーザは、予測に対するＡＲ成分の影響とＦＩＲ成分の影響とを調整することができる。
　ペアワイズで同定することの第三の利点は、第二の利点と密接に関係するが、各入力変数に対する（１）式のモデルの意味が明確であることである。ここで、先と同じ様に、出力変数を除く入力変数が１つの場合である下記（２）式の単純なＡＲＸモデルを考えてみる。
　　y(t)=a1×y(t-1)+ a2×y(t-2)+…+any×y(t-ny)
+b1×u(t-L)+ b2×u(t-L-1)+…+bnu×u(t-L-nu+1)　　（６）
　この時、上記（６）式をＡＲＸモデルとして同定すると、以下のｙに関する項（Ｙと定義する）とｕに関する項（Ｕと定義する）は、ＹとＵの意味は明確ではなく、ＹとＵとを個別に意味付けして説明することは困難である。
　　Y：=a1×y(t-1)+ a2×y(t-2)+…+any×y(t-ny)　　（７）
　　U：＝b1×u(t-L)+ b2×u(t-L-1)+…+bnu×u(t-L-nu+1)　　（８）
　一方、（７）式と（８）式とは、各々（５）式の形をしているので、ペアワイズで同定を行うと、ＹとＵとは各々ｙの予測値を意味することになるので、ＹとＵとは、各々別の説明変数（入力変数）に対する出力ｙの予測値であるという解釈をすることができる。
　このようにペアワイズの同定を行うと、（７）式と（８）式とから、（６）式の形の予測モデルを構成的（Constructive）に構築することは容易であり、例えば（７）式のＵと（８）式のＹとの平均（＝（Ｕ＋Ｙ）÷２）として構成することで、（６）式の形の予測モデルを得る事ができ、これは、再度ｙの予測値を意味していることは明らかである。
　上記のように、予測モデルの形（構造）が同じであっても同定の手順を変える事で、その意味づけは変化する。ペアワイズで同定することの利点は、意味付けのしやすいモデルを構築することができる点である。
　以上が多入力１出力ではなく１入力１出力のペアワイズで予測モデルを同定することの利点である。加えて、ペアワイズ予測モデルとして、（５）式のＦＩＲもしくはＡＲの形の回帰式の形で同定することで、ディジタル信号処理やシステム同定分野で従来から良く知られている周波数解析や安定解析などの様々な解析手法を、ペアワイズに直接流用できるため、その分野の知識があれば、ペアワイズでより詳細な解析や調整も可能になる。
　ペアワイズ予測モデルとして、（５）式の様なＡＲモデル／ＦＩＲモデルをさらに分解し、時間遅れを考慮した変数を各々一つの説明変数と見なした単回帰モデルに分解して予測モデルを構築することもできる。
　　y(t)=a×uk(t-L)+c　　（９）
　上記（９）式の様な単純な形にすると、（９）式に含まれるパラメータは、遅れ時間Ｌ、入力変数ｕｋに対する出力変数の倍率（ゲイン、比例係数）ａ、入力変数と出力変数の平均値の差（入力が０の場合の出力の値、バイアス）ｃの三つだけとなり、全てのパラメータの意味が明確になる。
　（５）式の形では、遅れ時間Ｌとバイアスｃの解釈は（９）式と同様であるが、ａｋ、ｋ＝１、２、…、ｎという回帰係数（パラメータ）の個々の意味までは明確でなく、システム同定やディジタル信号処理に関する知識を用いて、（５）式自身を解釈することはできるものの、個々のａｋの直感的な解釈は困難である。これに対し、（９）式の単回帰式の場合はａの解釈も容易であり、単純に、平均値を０に処理した場合の入力変数に対する出力変数の倍率を示している。これにより、全てのパラメータの意味を容易に解釈できる様になるため、（９）式を用いると、何等かのアルゴリズムを適用した場合のパラメータ同定結果が直感的な感覚とずれるような場合（例えば、異常データが混入している場合などにそのような事が起こり得る）、例えば予測モデルの構築者が直感的に手動でパラメータを調整することも容易になり、また、予測モデルの中身を他者に説明することも極めて容易になる。
　本実施形態のモジュラー型時系列データ予測装置において、例えば、出力変数はポンプ井流入量であり、入力変数として１か所の地上雨量であると仮定すると、（９）式のｃは地上雨量がゼロの場合の流入量を意味し、これは、汚水と雨水が同一管渠を通る合流式下水の場合は、汚水量に相当するものであるという解釈ができる。そして、（９）式のａは流入量から汚水量を引き去った場合に、地上雨量計による雨量に対する雨水流入量の倍率を示しており、地上雨量と流入量との相関が十分に高ければ、このような単純な倍率も意味を持つ。さらに、（９）式の遅れ時間Ｌは、地上雨量の観測点から雨水ポンプ井までの雨水の流れの遅れ時間を意味し、下水管に流出するまでの流出の遅れと下水管路を流れて雨水ポンプ井に到達するまでの流下時間の和を意味すると解釈することができる。
　このように、解釈や説明が極めて容易になることが（９）式の様な単純なモデルを用いることの利点であり、実際に（９）式の形だけを用いた最も簡単なモデルは、流入量予測の簡易モデルとして、特許文献１の中でも採用されている。
　上記のように直感的な解釈が可能になることにより、手動でのパラメータ調整が容易になる具体的な例について説明する。
　図３および図４は、幹線水位データと流入量データとの関係の一例を概略的に示す散布図である。
　図３および図４において、横軸は、幹線水位計１３１～１３Ｋのいずれかに相当するある箇所の幹線水位データ、縦軸は流入量計１１に相当する予測対象となる流入量データとした散布図の一例を示している。
　図３に示す直線は、最小２乗法によって、単回帰のゲインとバイアスを求めたものであり、この場合ゲインａ＝９．６０２４、バイアスｃ＝２３．３１４２となっている。これは、予測値の平均２乗誤差（ＭＳＥもしくはＲＭＳＥ）が最小になるという意味で、最適な回帰直線である。
　図３において、例えば幹線水位の値が－１．５であるとき、回帰直線に従うと流入量の値はおおよそ９程度と予測されることになる。しかしながら、図３の散布図によれば幹線水位の値が－１．５であるときの流入量は、９よりも大きい値である場合が圧倒的に多いことがわかる。
　上記のように、回帰直線による予測値と測定値との間にズレが生じる理由は、幹線水位が－２．５乃至－２付近のデータ数が圧倒的に多いため、この付近のデータにより適合する様に係数が推定されるためである。この回帰直線は、理論上は正しいものであるが、実際の流入量予測においては、必ずしも好ましくない場合が多い。なぜなら、流入量の予測値は、雨水排水ポンプの制御や運転支援に用いられるため、浸水リスクを回避するためには、流入量が多くなる場合をできる限り正確に予測したいという暗黙の要請があるためであり、流入量が少ない場合の予測精度はそれほど重要にならない場合が多いためである。
　しかし、このような暗黙の要請をアルゴリズムに組みこむ場合には、流入量が多い箇所の予測誤差に重みをつけるなどの処置が必要であり、その重みの調整なども含め、一般にかなりの労力を要する。
　図４では、最小２乗法と、最小２乗法の他の１０種類の様々な回帰手法とを、図３と同じデータに適用した結果得られた回帰直線Ｌ１－Ｌ１１の一例を示している。この中で、１番乃至１０番の方法は、Pythonと呼ばれるソフトウェアのパッケージに組み込まれている様々なアドバンストな回帰手法であるが、いずれの方法による回帰直線Ｌ１－Ｌ１０も、流入量が多いときに適切な予測値が得られるものではなかった。
　１１番目の方法は、幹線水位が－２．２乃至－１．５までの範囲のデータだけを抽出し、他のデータを使用せずに、通常の最小２乗法を適用した方法である。これは、本願発明者らが、幹線水位－２．５乃至－２付近のデータに適合するように回帰直線が求められることを理解した上で、その付近のデータを意図的に捨てることによって、流入量が多い場合のデータに適合するように意図的に調整した方法である。
　このように、アルゴリズムに頼った調整を行おうとすると、例えばデータの範囲を限定する、あるいは、データを間引くなどの様々な工夫を施さないと、本来望むパラメータが得られないことがある。しかし、１１番目の方法により得られる直線Ｌ１１は、２点を決めれば直線が唯一に定まるという極めて単純な原理を用いれば、人間が目視で直線Ｌ１１のようなラインを引くことにはほとんど労力を要しない。上記のように、単回帰という極めてわかりやすく、また図３や図４に示す散布図上に明確に可視化できるパラメータを用いることによって、手動での調整が極めて容易にできることがわかる。
　すなわち、単回帰モデルを用いると、倍率（ゲイン）とバイアスという容易に理解できるパラメータを用いて解釈できて、調整が容易になる。
　さらに、単回帰を用いると別の観点で解釈することも可能になる。すなわち、統計分野における基本的な統計量である、平均、標準偏差、および、相関、という三つのパラメータを用いて解釈することが可能になる。
　単回帰モデルの回帰係数（比例係数、傾き、倍率、ゲイン）ａは、相関係数ｒと次式の関係にあることが広く知られている。
　　a=r×(σy/σu)　　（１０）
　ここで、ｒ、σｙ、σｕは、各々、入力変数と出力変数の相関係数、出力変数の標準偏差、入力変数の標準偏差、である。また、回帰係数ａとバイアスｃは、単回帰モデルの予測２乗誤差が最小となる様に最小２乗法で同定した場合、以下の関係になることも広く知られている。
　　c=μy-a×μu　　（１１）
　ここで、μｙ、μｕは、各々、出力変数の平均と入力変数の平均、を表す。
　従って、（９）式において、入力変数と出力変数を各々、平均と標準偏差を用いて正規化しておけば、回帰係数ａと相関係数ｒの値は一致し、ｃ＝０となる。
　この関係が意味することは、単回帰モデルのゲインとバイアスでモデルを説明することと、各変数（入力変数と出力変数）の平均、標準偏差と、入力と出力との相関係数、のみでモデルを説明することは等価であるということである。従って、予め、変数毎に適切な平均と標準偏差の推定値を用いて各変数のデータを正規化しておけば、単回帰の入出力関係は、次式のように、相関係数のみで、全てを説明できることになる。
　　y´(t)=a×uk´(t-L)=r×uk´(t-L)　　（１２）
　ここで、y´=(y-μy)/σy、uk´=(uk-μuk)/σukで定義される正規化された出力変数と入力変数である。
　このように解釈すると、例えば、入出力の相関が無い場合、すなわちｒ＝０となる様な入出力関係を持つ場合は、その出力に対して対応する入力は全く予測能力を持たない事が（１１）式の関係から明確に理解することができる。また、現実のデータは、様々なノイズやアウトライア（外れ値）で汚染されている（コンタミされている）場合も多いが、このような時、平均、標準偏差、相関係数、の三つのパラメータを外れ値に対してロバストに推定するロバスト推定法を用いて推定することにより、外れ値に対してロバストな予測モデルを構築することもできる。もちろん、より複雑なモデルのパラメータ推定に対する様々なロバスト推定法も開発されているが、一般にロバスト推定した結果を解釈することは容易ではないのに対し、平均、標準偏差、相関係数という三つのパラメータに対するロバスト推定の場合は、何等かのロバスト推定を適用した推定結果の良否の判断や解釈も容易になるため、ノイズやアウトライアにコンタミされたデータに対する予測モデル構築を行いたい場合の解釈や調整も容易になる。
　このように、（９）式や（１０）～(１２)式は、容易に理解しやすいパラメータのみで構成されているが、（５）式と比較しても、極めて単純な構造をしており、実際の応用においては、その表現能力が十分でなく、十分な予測精度が得られない可能性があった。そこで、（９）式（あるいは（１０）～（１２）式）の形から、(５)式の形を合成することを検討する。（９）式のモデルを用いて(５)式の形を得るためには、(９)式の遅れ時間Ｌを可変にして足し合わせればよい。すなわち、まず、（９）式における遅れ時間Ｌ（（５）式のＬと区別するため以下ではＬ´とする）の最小の値を（５）式のＬと対応させてＬ´＝Ｌとし、Ｌ以上の遅れ時間について、Ｌ´＝Ｌ＋１、Ｌ´＝Ｌ＋２、…Ｌ´＝Ｌ＋ｎ－１としたｎ個の（９）式のモデルを加え合わせれば（５）式の形のモデルが得られる。
　すなわち、以下の（１３）式のｎ個の単回帰モデルを加え合わせれば、形式的に（５）式と等価な（１４）式のＦＩＲもしくはＡＲモデルの形のモデルが得られる。
　　y(t)=a1×uk(t-L)+c　　（１３＿１）
　　y(t)=a2×uk(t-L-1)+c　　（１３＿２）
　　y(t)=an×uk(t-L-n+1)+c　　（１３＿ｎ）
　　y(t)=1/n×(a1×uk(t-L)+a2×uk(t-L-1)+…+an×uk(t-L-n+1))+c　　（１４）
　異なるｎ個の遅れ時間のｎ個の単回帰モデルを作成した上で、その平均モデルを作成すれば、形式的には（１）式と同じＦＩＲあるいはＡＲモデルの形のモデルが得られる。この際、上記の様に単純に平均化処理を行うこともできるが、先の述べた単回帰のゲインとバイアスと相関係数との関係を用いて、以下の様に重み付き平均を行うこともできる。
　まず、予め、入出力変数を正規化して、（９）式を（１２）式の形で表しておく。すると、回帰係数は相関係数と一致するので、この回帰係数（の絶対値）で重みづけした重み付き平均値として（１４）式に類似した（１５）式のＦＩＲ／ＡＲモデルの形の式を得ることができる（以下では、正規化した変数と正規化していない変数との記号を区別せず、文脈に応じてｕ（ｔ）、ｙ（ｔ）は正規化した変数を表すこととする。）
　　y(t)=(a1×uk(t-L)+…+an×uk(t-L-n+1))/(|a1|+|a2|+…+|an|))　　（１５）
　このようにすると、相関が強くなる遅れ時間による予測に対して重みをつけたペアワイズ予測モデルの合成が可能になる。なお、回帰係数と相関係数の関係が（１２）式の様に陽に関係つけられない場合でも相関係数絶対値もしくは回帰係数の絶対値で重みづけ平均化処理を行うことは可能であるが、（１２）式の関係があることで、重みづけを行う事の意味がより明確で説得性の高いものとなる。また、（１５）式において、相関係数の絶対値ではなく相関係数の２乗によって重みづけを行っても良い。
　ペアワイズ予測モデルとして、上記の二つの合成法は、いずれの方法でも入出力間の線形の関係しか表現することができないが、入出力間に非線形の関係がある場合も多い。このような場合は、（５）式の代わりに、適当な非線形関数φ（・）を用いて、以下の（１６）式の様な非線形回帰の形として非線形変換すればよい。
　y(t)=a1×φ（uk(t-L)）+…+an×φ（uk(t-L-n+1)）+c　　（１６）
　これは、機械学習の分野で広く知られている非線形回帰のテクニックであり、パラメータａｋ、ｋ＝１、２、…、ｎに関する線形性さえ維持していれば、線形回帰の様々な手法を直接適用することができる。また、機械学習の分野でよく知られている様に、非線形関数φ（・）を直接指定しなくても、（１６）式を計算する際に必要となる計画行列と呼ばれる行列に対して、φ（・）を指定することと等価に変換できるカーネル関数（類似度関数）を直接指定しても良い。ただし、その場合には、直感的な解釈性が若干低下する可能性がある事には注意する必要がある。
　（１６）式の様な非線形回帰を用いる事で、ある入力変数と出力変数との間に非線形関係がある場合にも対応することが可能になる。なお、ペアワイズ予測モデルの同定では、入力変数毎に出力変数との回帰モデルを構築するので、ある入力変数と出力変数の関係は非線形であるが、別の入力変数と出力変数との関係は線形であるような場合には、入力変数毎に（５）式と（１６）式とを単に使い分けるだけでよく、これによって、一般的には説明性も向上することが期待される。以上の一連の作用が、本実施形態におけるペアワイズ予測モデル同定部４の作用である。
　次に予測モデル合成法定義部５では、ペアワイズ予測モデル同定部４で定義した予測モデルの合成法を定義する。
　最も簡単な予測モデルの合成法は、ペアワイズ予測モデルで各入力変数に対して出力変数を予測するモデルが構築されているため、その平均により予測モデルの合成を行う方法である。この場合は、予測モデル合成法定義部５では、ペアワイズ予測モデル同定部４のｐ個の予測出力（以下では、各入力変数ｕ１、ｕ２、…、ｕｐに対する予測出力をｙ１、ｙ２、…、ｙｐとする。）を平均化する次式を定義することになる。
　y(t)=mean(y1(t)、y2(t)、…、yp(t))=1/p×(y1(t)+y2(t)+…+yp(t))　　（１７）
　（１７）式の定義は、最も基本的な定義方法であるが、この定義を改良することにより、合成した予測出力に対する信頼性を向上させたり、意図的に予測に傾向（バイアス）を持たせた予測を行ったりすることが可能になる。これを以下に順に説明する。
　まず、（１７）式の様な単純な平均化処理（標本平均）を行うと、各ペアワイズ予測モデルの予測値が平等に扱われているため、予測精度の良いペアワイズ予測モデルの予測値と予測精度の悪いペアワイズ予測モデルの予測値が混合されて、合成した全体の予測精度が劣化してしまう可能性がある。また、実際の運用においては、ある入力変数の計測データの信頼性が低く、別の入力変数の計測データの信頼性が高いという場合も稀ではないが、信頼性の高い入力変数に対するペアワイズ予測モデルの予測値も信頼性の低い入力変数に対するペアワイズ予測モデルの予測値も平均化してしまうことで、合成した予測出力の予測精度が劣化してしまう。極端な場合、例えば、ある入力変数のセンサの故障や不具合などにより、当該入力変数の時系列データにアウトライア（外れ値、異常値）が多量に含まれる様になる場合などには、予測精度が劣化するだけでなくアウトライアに引きずられて合成した予測自身が無意味になる（破綻してしまう）可能性がある。
　このような状況に対応するために予測モデル合成法定義部５で、以下の様な複数のアプローチをすることができる。
　一つ目は、主に後者のアウトライアに対して合成した予測値自身が破綻しない事を重視するアプローチであり、（１７）式の標本平均処理に変えて、ロバスト推定を採用する方法である。
　（１７）式で用いる「標本平均」という処理は、ｙ１、ｙ２、…、ｙｐのｐ個の予測値の中から代表値を推定する方法の一つであり、このような代表値の推定を統計分野では位置母数の推定と呼ぶ。すなわち、「標本平均」は、位置母数の推定方法の一つであり、統計的には、標本平均は推定効率の観点では良い性質を持つことが知られているが、一方でアウトライアに対するロバスト性の面では最も脆弱（非ロバスト）であることも知られている。従って、外れ値に対するロバスト性を向上させるためには、（１７）式をロバストな位置母数推定で置き換える方が良い。
　ロバスト推定の分野では外れ値に対するロバスト性を評価するいくつかの指標が知られているが、最も直感的に理解しやすい指標としてブレークダウンポイントという指標がある。これは、統計的推定に用いるデータの中に何パーセントアウトライアが混入することを許容するかという指標であり、用いるデータのＸ％を、仮想的なアウトライアを想定して∞に置き換えた場合、推定量（平均などの位置母数を推定する場合は位置母数の推定値）が∞になる（＝破綻する）かならないかの境界（∞になる直前）のＸの値をブレークダウンポイントと呼ぶ。
　例えば、「標本平均」という位置母数の推定に対しては、ただ１点のデータを∞に置換するだけでその標本平均値も∞になるので、「標本平均」のブレークダウンポイントは０％である。ロバスト統計ではブレークダウンポイントの最大値は５０％であることが知られており、最大のブレークダウンポイントを持つ位置母数推定量として「中央値（メジアン）」が知られている。従って、予測値が外れ値に影響されないようにすることを最大の目的とする場合には、（１７）式の「平均」を「中央値」で置換して、予測モデル合成法定義部５の定義とすることもできる。
　一方、「中央値」はロバストではあるが、推定効率が悪い事が知られており、直感的にもｐ個の変数で予測している中の一つの予測値しか用いないため、精度の高い予測を行うことを重視する場合に最良の方法ではないことは容易に推測できる。そのため、推定効率の向上とロバスト性の向上を両立するための各種のロバスト推定方法が知られている。最も単純な方法は、「トリム平均（刈込平均）」と呼ばれる処理であり、ｐ個のデータの上位と下位とのα％を削除した上で平均をとる方法である。直感的に明らかな様にαを大きくするとロバスト性は向上し、その極限として中央値推定があり、αを小さくして０に近づけるとロバスト性が低下しその極限が通常の平均（標本平均）であることは明らかである。従って、このトリム平均を予測モデル合成法定義部５の定義とすることもできる。この場合、αの調整によってロバスト性を調整できるが、適切に調整する事自身が難しい場合も考えられるため、別のロバスト推定方法をとることもできる。
　このような方法の代表的な例として、ｐ個のｙ１、ｙ２、…、ｙｐの中からすべての組み合わせの（ｐ（ｐ－１）／２個）の、二つのｙｉとｙｊ（ｉ≠ｊ）との平均を計算した上で、その中央値を採用するホッジスレーマン推定量（ＨＬ推定）と呼ばれる位置母数推定を行うこともできる。この方法は、ブレークダウンポイントが約３０％でロバスト性が高いうえに推定効率も良いことが知られており、このＨＬ推定を予測モデル合成法定義部５の定義とすることもできる。
　他のロバスト推定法として、一般には多変量データに対して用いられる手法であるＭＣＤと呼ばれる手法が知られているがこれを適用することもできる。ＭＣＤは、ｐ個のデータの所定の割合（通常５０％～７５％）のデータを取り出し、そのすべての組み合わせの中で分散が最小になるデータを用いて推定を行う手法であり、この方法を用いて平均を推定することもできる。ＭＣＤもロバスト性が高く推定効率の良い方法として知られている。一方、ＨＬ推定やＭＣＤ推定は、ｐ個のデータの中から取り出すデータの数（ＨＬ推定場合２個、ＭＣＤの場合はｐ／２～３ｐ／４個程度）の全ての組み合わせに対する計算が必要となるため、ペアワイズモデルの入力変数の数ｐが増加すると処理時間が飛躍的に増加するため、予測モデル合成法定義部５の定義として用いても、実際に予測を行う際にリアルタイム性を確保できない可能性がある。
　このような場合に対応するため、「全ての組み合わせ」に対して処理を行うのではなく、リアルタイムで処理が可能な回数だけ繰り返し処理を行う様にブートストラップ法を用いて代用しても良い。すなわち、ｐ個のｙ１、ｙ２、…、ｙｐの中から、所定のｋ個のデータをランダムに繰り返し抽出し、所定の回数（Ｎ回）、その平均値を求める。繰り返し求めたＮ個の平均値の中から、ＨＬ推定と同じように、例えばＮ個の平均値の中央値を採用する。このようなブートストラップ法による平均値推定を、予測モデル合成法定義部５の定義とすることで、リアルタイム性を維持しながら、推定効率が良くロバスト性の高い予測値の合成を行うことができると考えられる。
　上記以外にも外れ値に重みを付けて推定を行うM推定などのロバスト推定を位置母数推定に適用した方法を予測モデル合成法定義部５の定義とすることもできる。
　二つ目は、先のアプローチが主に外れ値（外れた予測値）に対してロバストに予測値の代表値を推定することを目的としていたのに対し、実際に予測精度の良い予測値に重みを付けて合成した予測を行うことを目的としたアプローチである。基本的な考え方は、（１７）式の単純な標本平均に変えて、重み付き平均を行う方法であり、次式で表される式で合成予測出力を定義する。
　　y(t)=w1×y1(t)+w2×y2(t)+…+wp×yp(t)　　（１８）
　ここで、ｗｋ、ｋ＝１、２、…、ｐは、重みであり、ｗ１＋ｗ２＋…＋ｗｐ＝１となる制約を満たす。予測モデル合成法定義部５の定義では、この重みの設定法を定義する必要がある。
　これには、ペアワイズ予測モデル同定部４で、ペアワイズ予測モデルを同定した際の予測出力と実際の出力の間の（重）相関係数、あるいは、その２乗である決定係数を用いることができる。これにより、少なくともペアワイズの予測モデル同定時のデータに対して、予測精度の高い予測能力を持つモデルの予測出力を重視した重み付き平均で予測出力を合成することができる。
　また、ペアワイズの予測モデルを、例えば上述の（１５）式のように構築した場合には、各遅れ時間毎の説明変数と出力変数の間の相関係数が既に算出されているので、統合したペアワイズ予測モデル（（１４）式に相当）毎に、各相関係数の絶対値の平均値ｒｋ、ｋ＝１、２、３、…、ｐを求め、ｒｋに応じて、重みｗｋをｗｋ＝ｒｋ／（ｒ１＋ｒ２＋…＋ｒｐ）と定義することもできる。
　このように、ペアワイズ予測モデルの予測能力に応じて重み付き平均値の重みを決定して（１８）式により、予測モデルの合成法を定義する動作が、本実施形態における予測モデル合成法定義部５の動作の一例である。
　三つ目は、二つ目と同様に（１７）式の、各ペアワイズ予測モデルの重み付き平均値により、合成した予測出力を定義するが、この重みを、データを用いて同定（推定）することで決定するものである。この際、ペアワイズ予測モデルの同定時に用いたデータと同定したパラメータ値を用いると、その同定データに対する各ペアワイズ予測モデルの予測値ｙ１、ｙ２、…、ｙｐの時系列データが得られる。
　そして、（１７）式を、この予測値ｙ１、ｙ２、…、ｙｐを入力として、同定に用いた実際の出力データを出力ｙと見なすと、（１７）式自身が重みｗｋ、ｋ＝１、２、…、ｐを回帰係数とする重回帰モデルの形式になっている。従って、これらの予測時系列データと同定に用いた出力データとを用いて、重みｗｋ、ｋ＝１、２、…、ｐを重回帰の方法で同定することができる。ただし、重みは各々正でなければならないという不等式制約とｗ１＋ｗ２＋…＋ｗｐ＝１という等式制約を満たす必要があるため、通常の重回帰で用いる最小２乗法は適用できないが、混合線形推定法などの制約条件を考慮できる線形回帰の推定法を用いて、重みを推定することが可能になる。
　なお、この方法は、（１３）式から（１４）式を得る代わりに（１５）式を得る箇所にも同様の考え方を適用することができるので、単回帰モデル⇒ＦＩＲ／ＡＲモデル⇒多入力の伝達関数モデル、という３段階の段階的な学習（推定）によって予測モデルを合成することができる。
　本実施形態において、上記のようにして重みｗｋ、ｋ＝１、２、…、ｐを推定した（１７）式によって予測モデルの合成法を定義する方法が、予測モデル合成法定義部５の動作の他の例である。
　四つ目は、一つ目から三つ目とは異なり、合成した予測に傾向（バイアス）を持たせたい場合のアプローチであり、意図的に過大、あるいは、過小の予測を行う事を目的とする場合の予測出力の合成法である。
　上記の予測を行いたい状況は、現実の問題ではしばしば遭遇する。例えば、本実施形態で行う雨水の流入量予測の場合、この予測情報は雨水排水ポンプの起動や停止のタイミングを図るための支援情報として利用される場合が多い。このような場合、実際に流入する雨水流入量より過小な流入量の予測を行ってしまうことはリスク回避の観点から致命的になる事がある。
　このような場合、もちろん、正確な流入量を予測することが好ましい事は言うまでもないが、それが現実的に困難な場合は、多少多めの量を予測しておく方が安全である。
　また、例えば雨水排水ポンプの制御において流入量予測を利用する際には、ポンプの起動タイミングを図るために過大側の予測をすることが好ましいが、停止タイミングを図るためには過小側の予測をすることが好ましいため、過小側の予測を行いたい場合もある。
　また、例えば、水や電力の需要予測などの場合は、あまり過小に予測を行ってしまうと、水や電力の供給に支障をきたす可能性が考えられるため、若干過大側に予測を行っておく方が安全な場合もある。
　また、このように、過大側あるいは過小側の予測を行える様にしておくと、純粋に過大な予測や過小の予測を行いたいという場合以外にも、実際に予測システムを運用して長期評価した場合に、予測値が過小に出る傾向がある場合には過大側に再調整したり、逆に予測値が過大に出る場合には過小側に再調整したりして微調整を行うことも可能になる。
　このような動機に基づいて行う予測出力の合成法を以下に示す。これを行うためには、各ペアワイズ予測モデルの出力予測値ｙ１、ｙ２、…、ｙｐに対して、先に述べた位置母数の推定方法を定義する以外に、標準偏差や分散に対応する尺度母数の推定方法も同時に定義する。例えば、位置母数を平均とした場合には、尺度母数として、標準偏差を採用する様に定義する。また、位置母数としてロバスト推定法のメジアン（中央値）を用いた場合には、尺度母数としてＭＡＤ（中央値絶対偏差：Median Absolute Deviation）を採用する。さらに、ＨＬ推定、ＭＡＤ推定、ブートストラップ推定、Ｍ推定などを用いた場合には、各々の位置母数に対する尺度母数の推定方法を定義する。このようにして、ｙ１、ｙ２、…、ｙｐに対して適用する位置母数と尺度母数の推定法を予め定義しておく。これを以下では各々μとσと表す（通常μは平均値、σは標準偏差の意味で使われることが多いが、ここでは、必ずしも平均と標準偏差ではなく、ここで定義した位置母数と尺度母数を表すものとする）。
　次に、各ペアワイズ予測モデルの予測出力の集合、ｙ１、ｙ２、…、ｙｐに対して、その集合のμ±ｋσの値を合成した予測出力とするように、ｋの値と符号とを指定する。例えば、ｋ＝１としてμ＋σを指定すると、ｙ１、ｙ２、…、ｙｐの中の位置母数から尺度母数の１倍だけ過大な値を合成した予測出力とすることを意味する。逆にｋ＝－１とすると、ｙ１、ｙ２、…、ｙｐの中の位置母数から尺度母数の１倍だけ過小な値を合成した予測出力とすることを意味する。もし、位置母数と尺度母数として平均と標準偏差を用いた場合は、ｙ１、ｙ２、…、ｙｐが正規分布に従うと仮定すると、μ＋σは約６８％の値（μが５０％の値）、μ－σは約３２％の値を取り出すことに対応する。
　なお、平均や標準偏差を用いなくても、ｙ１、ｙ２、…、ｙｐの中のＫ％にあたるＫ％分位点を取り出して合成出力とするという方法も、位置母数と尺度母数を適切に定義すれば、μ±ｋσの形式で書けるので、このように設定した分位点と直接抽出するという方法で指定することも可能である。例えば、過大側の予測を行いたい場合に予測の異常値をある程度除外する事を想定して９０％分位点にあたる予測値を採用し過大予測を行う様に合成した予測出力を定義したり、１０％分位点にあたる予測値を採用し過小予測を行う様に合成した予測出力を定義したりすることもできる。
　予測モデル合成法定義部５の定義方法として四つの方法を説明したが、予測モデル合成法定義部５は、これらのいずれか、あるいは、その組み合わせによって予想モデルの合成法と定義することができる。なお、以上の一連の作用は、オフラインで過去のデータを用いて定期的もしくは非定期的に実行される。
　次に、オフラインで同定および定義した定義式を用いて、オンラインで予測を行う。　オンラインの予測は、時間の進行方向における所定の周期ＴＨ（＜＜ＴＬ（オフライン同定を定期的に行う場合の周期））で行われる。
　まず、オンライン予測用データ抽出部３２から周期ＴＨで、入力変数のデータを抽出する。本実施形態では、入力変数は、流量計（雨水ポンプ井流入量計）１１による雨水ポンプ井流入量と、幹線流量計１２による幹線流入量と、幹線水位計１３１～１３ＫによるＫ個の幹線水位計と、地上雨量計１４１～１４ＭによるＭ個の地上雨量と、Ｑ×Ｐメッシュ１５１１～１５ＱＰの各メッシュにおけるレーダ雨量と、を含む。
　次に、ペアワイズ出力変数予測部６では、ペアワイズ予測モデル同定部４で同定した各入力変数と出力変数のペアワイズ予測モデルを用いて、ｐ個の予測出力を計算する。これは、（５）式、（１４）式、（１５）式などの形で同定したパラメータ値が同定されたペアワイズ予測モデルにオンライン予測用データ抽出部３２で抽出したオンラインの入力変数データを入力することで直ちに計算できる。これがペアワイズ出力変数予測部６の作用である。
　次に、合成出力変数予測部７では、予測モデル合成法定義部５で定義した（１７）式や（１８）式に相当する合成法の定義式を用いて、ペアワイズ出力変数予測部６から出力される各入力変数に対するペアワイズの予測出力を入力することで合成した予測出力が計算される。
　次に、予測誤差評価部８は、合成出力変数予測部７で計算された予測出力の予測結果とペアワイズ出力変数予測部６のペアワイズの予測出力結果とを取得し、時間の進行方向における所定の周期ＴＨで時系列データとして保存する。ここで、保存された予測結果は、誤差評価を行うタイミングまで保持される。予測誤差の評価を行うタイミングは、予測モデルのユーザなどが外部から指示しても良いし、時間の進行方向における所定の周期ＴＭ（例えばＴＨ＜ＴＭ＜ＴＬ）で実施しても良い。いずれの方法であっても、予測誤差の評価を行うタイミングで、評価用データ抽出部３３で、予測誤差評価部８に保存された予測出力の時系列データと対応する時間の計測された出力変数、すなわち、本実施形態では雨水流入量の時系列データを抽出する。
　そして、予測誤差評価部８では、抽出された出力（雨水流入量）の時系列データとｐ個のペアワイズ予測モデルの予測出力との誤差、および、合成した予測出力の誤差を評価する。この際、誤差の評価は、通常時系列解析で行われるＭＳＥ（平方２乗誤差＝Ｌ２誤差）やその平方根であるＲＭＳＥなどであっても良いが、目的に応じて、誤差評価の評価基準を適宜設定しておいても良い。例えば、ピークの雨水流入量の予測精度が重要である場合はピーク流入量の差や誤差が最大になる場合の差（Ｌ∞誤差＝∞ノルム誤差）を誤差評価基準としても良いし、Nash-Sutcliffe係数などの評価基準を用いても良い。また、流入量が急激に増加する場合の遅れを改善する事を目的としたい場合には、予測出力と実績出力との位相差を計算し、位相差を小さくすることを誤差の評価基準としても良い。いずれにしても、予測誤差評価部８は、設定した誤差評価基準に基づいて、所定のタイミングもしくは周期ＴＭで、時系列データとｐ個のペアワイズ予測モデルの予測出力との誤差、および、合成した予測出力との誤差を評価する。
　次に、ペアワイズ予測モデル修正部９では、予測誤差評価部８で評価した予測誤差に基づいて、ペアワイズ予測モデルのパラメータの調整、あるいは、特定のペアワイズモデルの削除、もしくは、複数のペアワイズモデルの統合を行う。
　以下に、上記のような調整機能、削除機能、および、統合機能を付加したモジュラー型予測をより具体化した方法の例について説明する。
　一つ目のパラメータ調整について、ペアワイズ予測モデル修正部９は、予測誤差評価部８で評価した予測誤差に基づいて、著しく予測精度が悪い（予測誤差が所定の許容誤差より大きい）と判断される入出力ペアを抽出し、抽出された入出力ペアのペアワイズ予測モデルの再同定を行う。予測精度の悪さの判定は、（１）予測誤差に対するしきい値（許容誤差）を直接指定してしきい値を超過したペアワイズ予測モデルを予測精度の悪いモデル（あるいは入出力ペア）と判断する、（２）合成した予測出力の誤差とペアワイズ予測モデルの各予測誤差とを比較し、合成した予測誤差に対して相対的に著しく誤差が大きい（予め設定された許容範囲を超える）予測誤差を持つペアワイズ予測モデルを予測精度の悪いモデルと判断する、（３）ペアワイズ予測モデルの予測誤差（＋合成した予測誤差も含んでも良い）に対する位置母数μｅと尺度母数σｅを計算し、μｅ±ｋσｅ（ｋは設定するパラメータで２～３ぐらいで通常は設定する）から外れた予測誤差を持つペアワイズ予測モデルを予測精度の悪いモデルと判断する、などの方法で行うことができる。（１）－（３）のような判断によって、パラメータの再調整（再同定）が必要と判断された場合、ペアワイズ予測モデル修正部９は、所定の期間の直近のデータを用いてパラメータの再同定を行う。
　本実施形態において、例えば、Ｍ個の地上雨量データの中のいずれかの地上雨量データと雨水流入量との間のペアワイズ予測モデルの精度が悪化したと判断された場合、ペアワイズ予測モデル修正部９は、地上雨量と雨水流入量との間のペアワイズ予測モデルのパラメータのみを再同定する。このように部分的な同定による再調整を行う事ができる点は、モジュラー型の構成（構成的なアプローチ）の大きな利点である。
　二つ目の特定のペアワイズモデルの削除は、一つ目のパラメータの調整と類似の手続きで行うことができる。基本的な手順としては、ペアワイズ予測モデル修正部９は、一つ目のパラメータの調整で述べた方法と類似の方法で誤差評価を行い、誤差が著しく大きくなり精度が劣化したペアワイズ予測モデルを抽出し、抽出したペアワイズモデルを削除する。誤差評価の方法は一つ目の方法と同じであるが、ペアワイズ予測モデル修正部９は、（１）劣化の程度を判断するしきい値を一つ目のものよりも大きくして、一つ目の調整よりも明らかに予測精度が劣化しているものを削除対象とする、（２）まず一つ目のパラメータの再調整を行い、パラメータの再調整を行っても精度劣化の判断に用いた基準をクリアできない場合に削除対象とする、などの方法で削除すべきペアワイズ予測モデルを抽出する。
　本実施形態においては、例えば、幹線水位計などは水位計が水没することなどにより、水位を正確に測ることができなくなりアウトライア（異常値）が多量に混入してしまうことなどが考えられる。この場合、水没した水位計によるペアワイズ予測モデルは分離されることになる。このように、予測精度に悪影響を与えている部分を部分的に切り離すことで予測精度の改善を図れることもモジュラー型の構成（構成的なアプローチ）の大きな利点である。
　３つ目のペアワイズ予測モデルの統合は、以下の様に実施される。まず、ペアワイズ予測モデル修正部９は、ペアワイズ予測モデルの予測誤差を相互に比較評価し、その予測誤差同士の差の絶対値が所定のしきい値以下となる相互に類似の予測精度を持つペアワイズ予測モデルを、相互に類似する予測モデルとして抽出しグループ化する。
　次に、ペアワイズ予測モデル修正部９は、このようにしてグループ化されたペアワイズ予測モデルの持つパラメータの値をグループ内で相互に比較する。例えば、係数パラメータａｋ、ｋ＝１、２、…、ｎをベクトル化してＡ＝［ａ１、ａ２、…、ａｎ］とし、相互に比較する予測モデルの係数パラメータベクトル同士の類似度を例えば内積などにより算出し、所定のしきい値により類似していると判断されたペアワイズ予測モデルを統合の対象候補とする。
　本実施形態のモジュラー型の構成的なアプロ―チでは、パラメータの可同定性（一意決定可能性）が基本的に担保されているので、予測結果が類似しておりパラメータ値も類似していれば、その入力変数はほぼ同じ値を持つはずであるが、念のため、ペアワイズ予測モデル修正部９は、候補として抽出されたペアワイズ予測モデルの入力変数同士の相関係数を求めて、所定のしきい値（例えば、相関係数０．９５）以上である場合に最終的に統合すべきペアワイズ予測モデルと判断する。そして、ペアワイズ予測モデル修正部９は、最終的に統合すべきと判断された複数のペアワイズ予測モデルの入力変数の平均をとった変数を、新たな一つの合成された入力変数として定義する。
　そして、ペアワイズ予測モデル修正部９は、（１）先に比較した係数ベクトルの平均値を採用する、あるいは、（２）新たに合成した平均化された入力変数と出力変数に対して再同定を行って推定する、のいずれかの方法を用いて、合成された入力変数に対するペアワイズ予測モデルの係数を決定する。
　上記のような統合は、本実施形態のモジュラー型時系列データ予測装置においては、例えば次の様なケースにおいて生じると考えられる。入力変数の中にＱ×Ｐメッシュ１５１１～１５ＱＰの各々のレーダ雨量のデータが含まれているケースにおいて、レーダ雨量のメッシュのサイズは、例えばXRAINと呼ばれるＸバンドの気象レーダの場合約２５０ｍ×２５０ｍであるため、隣接するメッシュの雨量データはほとんど同じ値をとる場合が多い。このような場合、相互に隣接するメッシュの雨量データは統合の対象となる可能性が高く、Ｑ×Ｐ個の入力変数が、少数の互いに独立した合成入力変数として統合されることになる可能性が高い。このような統合が起こると、レーダ雨量のメッシュ状のデータが、雨水流入量に与える影響の類似性に応じて統合されるため、入力変数の数を減らせるだけでなく、説明性が向上する事が期待できる。
　また、このような場合は、先に述べた多重共線性の問題が生じているケースであるが、これに対し、例えば、Lassoと呼ばれる方法などを用いて入力変数を自動的に選択すると、いくつかのメッシュの雨量データが代表的な入力変数として選択され、選択された入力変数に対応するメッシュの雨量と極めて類似するメッシュの雨量データは無視されることになる。雨量データが極めて類似しているので、このように代表メッシュの雨量を入力変数として選択しても通常は問題がない。しかしながら、近年頻発している局所的な豪雨が見られる場合、例えば、たまたま選択したメッシュの雨量はあまり多くないのに、入力変数として選択されなかった近くのメッシュの雨量が非常に大きい様な局所的な降雨が生じた場合、あるメッシュを代表メッシュとして選択してしまうと、流入量が増加しないという誤った予測を行う可能性がある。
　一方、ここで述べた入力変数の統合を行うと、類似の影響を与えるメッシュの雨量が平均化され、これは統合されたメッシュの平均雨量という意味を持つため、このようなケースでも妥当な予測結果を与える事が期待できると同時に、合理的な説明性も向上する。このように統合すべき入力変数をまとめて解釈し、説明しやすくするようにできる事もモジュラー型の構成（構成的なアプローチ）の大きな利点である。
　最後に、合成出力変数予測部７で算出（合成）された予測出力は、出力予測結果観測部１０に送られ、例えばいわゆるトレンドグラフと呼ばれる時系列データとしてユーザであるプラントの管理者や運転員に提示される。この際、リアルタタイム監視している監視の現在までの値は、実際の出力値である雨水流入量とその予測値を同時に表示し、将来の予測値は予測出力のみを表示する様にしておくことが好ましい。
　また、合成した予測出力（予測雨水流入量）だけでなく、ペアワイズ予測モデルのｐ個の予測出力を用いて、その位置母数μｙと尺度母数σｙを計算し、μｙ±Ｋσｙ（Ｋは設定パラメータで２～３程度で設定する）の範囲を同時表示したり、あるいは、ペアワイズ予測モデルの予測出力の中から最大値（あるいはロバスト性を考えて最大値に近い９５％点にもっとも近い予測値）や最小値（あるいはロバスト性を考えて最小値に近い５％点にもっとも近い予測値）を直接予測出力の範囲として同時に表示したりしても良い。
　本実施形態のモジュラー型時系列データ予測装置によれば、一般にはブラックボックスモデルと呼ばれるモデルに対して、その内部パラメータに対して物理法則に矛盾しない合理的な解釈が可能になる様な仕組みを導入することで、ホワイトボックスモデルに近い、モデルの解釈と合理的な説明性を向上させて、部分的な調整、削除、統合などのモデルの維持管理持管理（メンテナンス）の煩雑さを大幅に低減できるという効果が得られる点である。
　なお、このような効果を得るための根本的なアイデアは、本来、入出力データしか与えられていない内部構造が完全に未知である状況下で、出力に関連する可能性のある要素を要素毎に出力と直接関係づけ（ペアワイズ予測モデルに対応）、各要素と出力の関係を、最終的な入出力関係の内部構造として定義される様に、構成的に予測モデルを構築していることである。このようなアイデアを用いることで、本実施形態のモジュラー型時系列データ予測装置では、本来知ることのできない内部構造に立ち入ることなく、疑似的にホワイトボックスモデルと同様の説明性と調整の容易さを実現している。
　また、このモジュラー型時系列予測モデルの構築方法は、統計的機械学習分野でよく知られている。先に述べたバギングやブースティングなどのアンサンブル学習と呼ばれる手法とも類似している。バギングやブースティングは、弱学習器（弱学習モデル）と呼ばれる、簡単に構成できるが、精度が必ずしも高くないモデルをたくさん集めることで、精度の高い強学習器（強学習モデル）と呼ばれる予測モデルを構築する方法であり、バギングは弱学習器に対して多数決などの処理で強学習を行う方法、ブースティングは、弱学習器の結果を用いて強学習器を学習させる方法である。ただし、バギングやブースティングでは、学習(同定)に利用データを変えることでたくさんの弱学習器を生成しているが、本実施形態では、同定データを変えるのではなく、説明変数を変えることで、弱学習器に相当するペアワイズ予測モデルを同定している点で異なっている。なお、本実施形態のモジュラー型時系列データ予測装置において、合成出力の予測方法を定義する際に、ロバスト位置母数推定や予め指定した重みを用いた重み付き平均を用いるような方法はバギングの考え方に近く、重み付き平均の重みを再同定によって求める方法はブースティングの考え方に近い。
　以下、本実施形態のモジュラー型の時系列予測装置におけるペアワイズ予測モデルの複数の実施例について説明する。
（第１実施例）
　図５は、一実施形態のモジュラー型の時系列予測装置におけるペアワイズ予測モデルの第１実施例について説明するための図である。
　図５では、（１）式のＦＩＲモデルやＡＲモデルをペアワイズ予測モデルとした例を示している。
　本実施例のペアワイズ予測モデルＭＡは、回帰ブロックＡ１－ＡＮと、予測出力の合成ブロックＡ５と、を備えている。
　例えば回帰ブロックＡ１－Ａ３のそれぞれは、入力変数に対して重み付き移動平均成分（ＦＩＲ成分／ＭＡ成分）を予測値として出力する。例えば回帰ブロックＡＮは、入力変数（過去の出力変数）に対して自己回帰成分（ＡＲ成分）を予測値として出力する。
　予測出力の合成ブロックＡ５は、回帰ブロックＡ１－ＡＮから出力された予測値を取得し、予測値を合成した合成予測値（Y(t+1)=F(Y1(t+1),Y2(t+1),…,Yn(t+1))）を演算して出力する。
　例えば、通常のＡＲＸモデルなどによってパラメータを同定すると自己回帰成分（ＡＲ成分）の重みが大きくなり、ＦＩＲ成分に相当する入力変数の値が急激に変化した場合に、予測出力が計測出力に追従できなくなり予測に位相遅れが生じる事がある。しかし、本実施例の構成をとることにより、ＡＲ成分の影響を調整することが可能になり、位相遅れを改善して計測出力の値が変化する前に出力の変化を予測することが可能になる。
（第２実施例）
　図６は、一実施形態のモジュラー型の時系列予測装置におけるペアワイズ予測モデルの第２実施例について説明するための図である。
　図６は、合成した予測に傾向（バイアス）を持たせたい場合に用いられる予測モデルの一例として、意図的に過大の予測および過小の予測を行う場合の予測モデルを示している。
　本実施例のペアワイズ予測モデルＭＢは、回帰ブロックＢ１－ＢＮと、予測出力の合成ブロックＢ５と、を備え、合成ブロックＢ５が、過大推定予測値（Y+(t+1)=F+(Y1(t+1),Y2(t+1),…,Yn(t+1))）と、通常推定予測値（Y(t+1)=F(Y1(t+1),Y2(t+1),…,Yn(t+1))）と、過小推定予測値（Y－(t+1)=F－(Y1(t+1),Y2(t+1),…,Yn(t+1))）と、を出力する点において上述の第１実施例のペアワイズ予測モデルＭＡと異なっている。
　本実施例では、ペアワイズ予測モデルの位置母数と尺度母数とを用いて、通常の予測に加えて過大な予測および過小な予測を行うことにより、安全側の予測（リスクを回避する側の予測）などのニーズ（要求）に応じて、予測結果を微調整することができる。
　図７および図８は、第２実施例のペアワイズ予測モデルの効果を説明するための図である。
　図７は、回帰手法の一種であるL1正則化を導入したLasso回帰による雨水流入予測値と測定値との一例を示している。
　図８は、第２実施例の予測モデルによる雨水流入予測値と測定値との一例を示している。
　図７に示す雨水流入予測値と測定値とを比較すると、予測誤差は小さいが、若干推定値が測定値よりも小さい結果となっており、過小推定になっている事がわかる。
　一方、図８に示す本実施例の手法で得られた通常の推定（平均的な推定）、過大推定、および、過少推定を行った結果と測定値とを比較すると、通常の推定値（予測流入量）では予想誤差は小さいく、ピーク付近での予測精度は良いものの流入量が増加するタイミングでの予測値が測定値を下回っている。
　このような場合、この予測値を制御に用いる場合には、もう少し流入量急増時の予測を過大に評価しておきたい場合がある。このとき、過大推定予測を行うと、全体的にＭＳＥ（平均２乗誤差）での予測精度は劣化するが、流入量急増時に急増する可能性を示唆する予測が可能になることがわかる。また、通常の予測値よりも過小に評価をしておきたい場合には、過小推定予測を行うと、予測値が小さくなるように調整することができる。
　このような過大推定および過小推定により予測値を調整可能とすると、長期の運用の結果、通常の予測値が実測値に対して過小推定、あるいは、過大推定になっている様な傾向（バイアス）が見られる場合、通常の予測値を過大推定値、あるいは過小推定値に置き換えることで、バイアスを補正してより精度の高い予測を行うことができる。
（第３実施例）
　図９は、一実施形態のモジュラー型の時系列予測装置におけるペアワイズ予測モデルの第３実施例について説明するための図である。
　本実施例では、ペアワイズ予測モデルとして、（５）式のようなＡＲモデル／ＦＩＲモデルをさらに分解し、時間遅れを考慮した変数を各々一つの説明変数と見なした単回帰モデルに分解して予測モデルを構築している。
　本実施例のペアワイズ予測モデルＭＣは、複数の時間遅れブロックと、複数の単回帰ブロックと、相関重み付き平均ブロックと、を含む予測出力ブロックＣ１－ＣＮを入力変数１－Ｎ毎に含み、当該複数の予測出力ブロックＣ１－ＣＮから出力された予測値を合成した最終合成予測値を出力する予測出力値の合成ブロックＣ（Ｎ＋１）を更に備える。
　予測出力ブロックＣ１は、複数（ｍ）の時間遅れブロック１１１－１１ｍと、複数（ｍ）の単回帰ブロック１２１－１２ｍと、相関重み付き平均ブロック１３と、を含む。
　複数の時間遅れブロック１１１－１１ｍは、入力変数１（例えば時刻ｔにおける値）に対してそれぞれに割り当てられた遅れ時間Ｌ１１－Ｌ１ｍ１だけ遅れた変数（例えば時刻ｔ－Ｌ１１、…、ｔ－Ｌ１ｍ１における値）を出力する。
　複数の単回帰ブロック１２１－１２ｍは、単回帰係数を用いて、時間遅れブロック１１１－１１ｍから出力された変数に基づく予測値を演算する出力する。
　相関重み付き平均ブロックは、複数の単回帰ブロック１２１－１２ｍから出力された複数の予測値の相関重み付き平均値を演算して出力する。
　予測出力ブロックＣ２－ＣＮは、上記予測出力ブロックＣ１と同様の構成であるため、個々の説明は省略する。
　予測出力値の合成ブロックＣ（Ｎ＋１）は、予測出力ブロックＣ１－ＣＮから出力された予測値を合成（例えば、相関重み付き平均）して、最終推定予測値（Y(t+1)=F(Y1(t+1),Y2(t+1),…,Yn(t+1))）を出力する。
　本実施例のペアワイズ予測モデルは、遅れ、ゲイン、バイアスの３項でモデルを説明することが可能であり、パラメータの微調整が容易である。また、単回帰係数と相関係数とには陽な関係があるため、各入力変数と各出力変数との相関係数を求めるだけで全体の予測モデルの構築ができ、予測モデルの解釈が容易である。
　更に、本実施例のペアワイズ予測モデルでは、第１実施例のペアワイズ予測モデルと同様の効果を得ることができる。
　図１０および図１１は、第３実施例のペアワイズ予測モデルの効果の一例を説明するための図である。
　図１０および図１１では、第３実施例のペアワイズ予測モデルに対し１０個の入力変数を入力して雨水流入予測を行った例である。ここでは、各入力変数に対する遅れ時間を一つだけ考えた最も単純な場合の予測結果の例を示している。
　入力変数の一つである降雨強度（地上雨量）の遅れ時間を推定すると３３分という推定結果が得られた。図１１は、上記推定結果を用いた雨水流入量の予測結果の一例を示している。図１１の予測結果によれば、全体的には、ある程度の精度で予測ができているが、流入量が急増する部分での予測の遅れ（位相遅れ）が大きい事がわかる。
　この時、位相遅れが大きい理由として、遅れ時間３３分というものが大きすぎる可能性が推測できる。この推測に基づき、図１０に、遅れ時間を１５分と調整した場合の雨水流入量の予測結果の一例を示す。図１０に示す予測結果によれば、遅れ時間を調整することにより、全体的な予測精度には大きな影響を与えることなく、流入量急増時の位相遅れを大幅に改善できていることがわかる。
　このように、本実施例のペアワイズ予測モデルによれば、容易に解釈可能なパラメータのみで予測モデルを構築可能であり、説明性が向上するだけでなく予測モデルの微調整を容易に行うことができる。
　また、単回帰係数と相関係数とは、入出力データを各々平均と標準偏差とで正規化すれば、一致することが知られており、これを逆に非正規化すれば、回帰係数とバイアスとに換算することができる。したがって、各入出力データの平均、標準偏差、および、入出力の相関のみを計算すれば、遅れ時間を指定するだけで、全ての予測モデルのパラメータを同定することができる。
　さらに、合成出力の計算を相関係数の重み付き平均で計算する方法とすることにより、予測モデルを構築するために必要な計算は、遅れ時間の推定（指定）と、入出力変数の平均および標準偏差と、入出力変数の相関と、のみであり、予測モデルに含まれるパラメータの全てを、この４種類の推定値（遅れ時間、平均、標準偏差、相関係数）だけを用いて説明することが可能になる。
（効果）
　上記のように、本実施形態によれば、入力変数と出力変数との相関関係が変化したり劣化したり、入力変数の時系列データに多量のアウトライアが含まれたりするなどの理由で予測精度が十分でない場合に、全体の予測モデルを再同定することなく、合成した予測出力値に影響を与えるペアワイズの予測モデルの部分的な調整や、悪影響を与える部分の削除を容易に行う事ができると同時に、これにより、予測精度劣化の理由を合理的に説明することが可能になる。
　また、入力変数間の多重共線性が極めて強くほぼ同じ説明能力を持つ入力変数が含まれる場合にも、パラメータの可同定(一意決定可能性)を維持することができ、物理法則に矛盾しない合理的な説明が可能な予測モデルを構築できると同時に、ほぼ同じ説明能力を持つ入力変数を統合することで、予測精度を維持しながら、説明性を向上させることができる。
　すなわち、本実施形態によれば、予測結果を合理的に説明するとともに、容易に予測結果を調整可能とするモジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラムを提供することができる。
　本発明のいくつかの実施形態を説明したが、これらの実施形態は、例として提示したものであり、発明の範囲を限定することは意図していない。これら新規な実施形態は、その他の様々な形態で実施されることが可能であり、発明の要旨を逸脱しない範囲で、種々の省略、置き換え、変更を行うことができる。これら実施形態やその変形は、発明の範囲や要旨に含まれるとともに、特許請求の範囲に記載された発明とその均等の範囲に含まれる。

Claims (12)

  1. 　複数のプロセス変数を所定の周期で計測する複数のプロセスセンサを有するシステム又はプロセスに適用される装置であって、
    　複数の前記プロセス変数の時系列データを所定の周期で収集し保存するとともに、複数の前記プロセス変数の中から予測対象となる少なくとも一つの出力変数を選択する出力変数データ選択部と、複数の前記プロセス変数の中から複数の入力変数の候補を選択する入力変数データ選択部と、を含むデータ収集保存部と、
    　前記時系列データから抽出された前記出力変数と複数の前記入力変数との同定用データを用いて、１入力１出力のペア毎にペアワイズ予測モデルのパラメータを同定して複数の前記ペアワイズ予測モデルを定義するペアワイズ予測モデル同定部と、
    　複数の前記ペアワイズ予測モデルから出力されるペアワイズ予測値の合成法を定義する予測モデル合成法定義部と、
    　時間の進行方向に所定の周期又はリアルタイムで前記時系列データから抽出された複数の前記入力変数の予測用データを、複数の前記ペアワイズ予測モデルに入力して、複数の前記入力変数のそれぞれに対応する前記ペアワイズ予測値を演算するペアワイズ出力変数予測部と、
    　前記合成法により複数の前記ペアワイズ予測値を合成して前記出力変数の予測値を演算する合成出力変数予測部と、
    　を有するモジュラー型時系列データ予測装置。
  2. 　前記ペアワイズ予測モデル同定部は、前記出力変数以外の複数の前記入力変数の同定用データに有限インパルス応答モデルを適用し、前記出力変数に自己回帰モデルを適用する、請求項１記載のモジュラー型時系列データ予測装置。
  3. 　前記ペアワイズ予測モデル同定部は、複数の前記入力変数の各々ついて１又は複数の遅れ時間を組み込んだ単回帰モデルを用いて、１入力１出力の伝達関数モデルを構築し、前記単回帰モデルを合成したモデルを前記ペアワイズ予測モデルとする、請求項１記載のモジュラー型時系列データ予測装置。
  4. 　前記ペアワイズ予測モデル同定部は、複数の前記入力変数各々の同定用データに対して非線形変換を行い、前記ペアワイズ予測モデルを構築する、請求項１記載のモジュラー型時系列データ予測装置。
  5. 　前記予測モデル合成法定義部は、複数の前記ペアワイズ予測値の重み付き平均を演算することを前記合成法とし、前記合成法における重みを、前記ペアワイズ予測モデル同定部によって同定された複数の前記ペアワイズ予測モデル各々の予測精度を表す指標に基づいて決定する、請求項１記載のモジュラー型時系列データ予測装置。
  6. 　前記予測モデル合成法定義部は、複数の前記ペアワイズ予測値の重み付き平均を演算することを前記合成法とし、前記ペアワイズ予測モデル同定部で用いた同定用データと同じデータを用いて、前記合成法における重みの値を同定して決定する、請求項１記載のモジュラー型時系列データ予測装置。
  7. 　前記ペアワイズ出力変数予測部から出力される複数の前記ペアワイズ予測値と、前記合成出力変数予測部から出力される前記出力変数の予測値とに対して、予め指定した所定の周期毎あるいは所定期間で、前記時系列データから抽出した実績値との誤差を評価する予測誤差評価部と、
    　複数の前記ペアワイズ予測値の誤差と、前記出力変数の予測値の誤差とを比較することにより、調整すべき前記ペアワイズ予測モデル、統合すべき前記ペアワイズ予測モデル、および、分離すべき前記ペアワイズ予測モデルを決定して、前記ペアワイズ出力変数予測部で用いられる複数の前記ペアワイズ予測モデルを修正するペアワイズ予測モデル修正部と、
    　を更に有する請求項１記載のモジュラー型時系列データ予測装置。
  8. 　前記ペアワイズ予測モデル修正部は、前記予測誤差評価部によって評価された複数の誤差の少なくとも一つが所定の許容誤差を超過したときに、当該誤差に対応する前記入力変数の候補に関する前記ペアワイズ予測モデルの再同定を行い、前記ペアワイズ出力変数予測部で用いられる複数の前記ペアワイズ予測モデルを修正する請求項７記載のモジュラー型時系列データ予測装置。
  9. 　前記ペアワイズ予測モデル修正部は、前記予測誤差評価部によって評価された複数の前記誤差の少なくとも一つが所定の許容誤差を超過したときに、当該誤差に対応する前記入力変数に関する前記ペアワイズ予測モデルを削除して、前記ペアワイズ出力変数予測部で用いられる複数の前記ペアワイズ予測モデルを修正する請求項７記載のモジュラー型時系列データ予測装置。
  10. 　前記ペアワイズ予測モデル修正部は、前記予測誤差評価部によって評価された複数の前記誤差を比較し、複数の前記誤差の差が所定の値以下であるときに、当該誤差に対応する前記入力変数の和あるいは平均値を一つの新たな入力変数として定義し、前記ペアワイズ出力変数予測部で当該新たな入力変数に関する前記ペアワイズ予測モデルが用いられるように、複数の前記ペアワイズ予測モデルを修正する請求項７記載のモジュラー型時系列データ予測装置。
  11. 　複数のプロセス変数を所定の周期で計測する複数のプロセスセンサを有するシステム又はプロセスに適用される方法であって、
    　複数の前記プロセス変数の時系列データを所定の周期で収集して保存し、
    　複数の前記プロセス変数の中から予測対象となる少なくとも一つの出力変数を選択し、　複数の前記プロセス変数の中から複数の入力変数の候補を選択し、
    　前記時系列データから抽出された前記出力変数と複数の前記入力変数との同定用データを用いて、１入力１出力のペア毎にペアワイズ予測モデルのパラメータを同定して複数の前記ペアワイズ予測モデルを定義し、
    　複数の前記ペアワイズ予測モデルから出力されるペアワイズ予測値の合成法を定義し、　時間の進行方向に所定の周期又はリアルタイムで前記時系列データから抽出された複数の前記入力変数の予測用データを、複数の前記ペアワイズ予測モデルに入力して、複数の前記入力変数のそれぞれに対応する前記ペアワイズ予測値を演算し、
    　前記合成法により複数の前記ペアワイズ予測値を合成して前記出力変数の予測値を演算する、
    　モジュラー型時系列データ予測方法。
  12. 　コンピュータに、請求項１１に記載の方法を実行させるモジュラー型時系列データ予測プログラム。
PCT/JP2022/018118 2021-05-10 2022-04-19 モジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラム WO2022239609A1 (ja)

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2021-079844 2021-05-10
JP2021079844A JP2022173863A (ja) 2021-05-10 2021-05-10 モジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラム

Publications (1)

Publication Number Publication Date
WO2022239609A1 true WO2022239609A1 (ja) 2022-11-17

Family

ID=84029588

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2022/018118 WO2022239609A1 (ja) 2021-05-10 2022-04-19 モジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラム

Country Status (2)

Country Link
JP (1) JP2022173863A (ja)
WO (1) WO2022239609A1 (ja)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116090388A (zh) * 2022-12-21 2023-05-09 海光信息技术股份有限公司 芯片内部电压预测模型生成方法、预测方法及相关装置
CN117036112A (zh) * 2023-10-09 2023-11-10 石家庄坤垚科技有限公司 一种土地规划用的地理信息系统及方法
CN118014719A (zh) * 2024-04-08 2024-05-10 南京启尚数字科技有限公司 一种基于线性回归模型的企业信用智能分析方法和系统

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
KR102640983B1 (ko) * 2022-12-21 2024-02-23 재단법인차세대융합기술연구원 극단치와 증감 추세를 반영하여 시계열 데이터를 기호화하는 분석 서버 및 그것의 데이터 분석 방법

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2004062440A (ja) * 2002-07-26 2004-02-26 Toshiba Corp 予測モデルシステム
JP2005222444A (ja) * 2004-02-09 2005-08-18 Toshiba Corp 統計的予測値演算方法および装置
JP2012515902A (ja) * 2009-01-21 2012-07-12 ボッカー，セバスチアン マススペクトロメトリーにより特に未知物質を同定するための方法
JP2016507113A (ja) * 2013-02-05 2016-03-07 ヨコガワ・コーポレーション・オブ・アメリカ 製品またはプロセス流の特性を判定するためのシステム、方法、および装置
JP6261960B2 (ja) * 2013-11-13 2018-01-17 株式会社東芝 雨水排水ポンプ制御装置、雨水排水システム、および雨水排水ポンプ制御プログラム
JP2018116545A (ja) * 2017-01-19 2018-07-26 オムロン株式会社 予測モデル作成装置、生産設備監視システム、及び生産設備監視方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2004062440A (ja) * 2002-07-26 2004-02-26 Toshiba Corp 予測モデルシステム
JP2005222444A (ja) * 2004-02-09 2005-08-18 Toshiba Corp 統計的予測値演算方法および装置
JP2012515902A (ja) * 2009-01-21 2012-07-12 ボッカー，セバスチアン マススペクトロメトリーにより特に未知物質を同定するための方法
JP2016507113A (ja) * 2013-02-05 2016-03-07 ヨコガワ・コーポレーション・オブ・アメリカ 製品またはプロセス流の特性を判定するためのシステム、方法、および装置
JP6261960B2 (ja) * 2013-11-13 2018-01-17 株式会社東芝 雨水排水ポンプ制御装置、雨水排水システム、および雨水排水ポンプ制御プログラム
JP2018116545A (ja) * 2017-01-19 2018-07-26 オムロン株式会社 予測モデル作成装置、生産設備監視システム、及び生産設備監視方法

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116090388A (zh) * 2022-12-21 2023-05-09 海光信息技术股份有限公司 芯片内部电压预测模型生成方法、预测方法及相关装置
CN116090388B (zh) * 2022-12-21 2024-05-17 海光信息技术股份有限公司 芯片内部电压预测模型生成方法、预测方法及相关装置
CN117036112A (zh) * 2023-10-09 2023-11-10 石家庄坤垚科技有限公司 一种土地规划用的地理信息系统及方法
CN117036112B (zh) * 2023-10-09 2023-12-22 石家庄坤垚科技有限公司 一种土地规划用的地理信息系统及方法
CN118014719A (zh) * 2024-04-08 2024-05-10 南京启尚数字科技有限公司 一种基于线性回归模型的企业信用智能分析方法和系统

Also Published As

Publication number Publication date
JP2022173863A (ja) 2022-11-22

Similar Documents

Publication Publication Date Title
WO2022239609A1 (ja) モジュラー型時系列データ予測装置、モジュラー型時系列データ予測方法、および、プログラム
US11238356B2 (en) Method of predicting streamflow data
CN111222698B (zh) 面向物联网的基于长短时记忆网络的积水水位预测方法
Mehdizadeh et al. Hybrid artificial intelligence-time series models for monthly streamflow modeling
Jung et al. Current status and future advances for wind speed and power forecasting
Belayneh et al. Standard precipitation index drought forecasting using neural networks, wavelet neural networks, and support vector regression
Sahoo et al. Application of support vector regression for modeling low flow time series
US20150206427A1 (en) Prediction of local and network-wide impact of non-recurrent events in transportation networks
KR102446493B1 (ko) 확률분포에 기반한 가뭄 예측 방법 및 이를 위한 장치
Kim et al. Comparative analysis of long short-term memory and storage function model for flood water level forecasting of Bokha stream in NamHan River, Korea
JP2007205001A (ja) 流量予測装置
CN103810532B (zh) 优化城市排水系统运行状况的方法
US20210088369A1 (en) Blockage detection using machine learning
JP4102131B2 (ja) 予測モデルシステム
CN104090974A (zh) 展延水库后续来水的动态数据挖掘方法及系统
Kasiviswanathan et al. Flood frequency analysis using multi-objective optimization based interval estimation approach
Wu et al. Multiple-clustering ARMAX-based predictor and its application to freeway traffic flow prediction
Ouyang Input optimization of ANFIS typhoon inundation forecast models using a Multi-Objective Genetic Algorithm
Tounsi et al. On the use of machine learning to account for reservoir management rules and predict streamflow
Zhou et al. Prediction and early warning method of inundation process at waterlogging points based on Bayesian model average and data-driven
Sahoo et al. Performance comparison of LS-SVM and ELM-based models for precipitation prediction in Barak valley: A case study
Turki et al. A Markova-Chain Approach to Model Vehicles Traffic Behavior
JP4423607B2 (ja) 降雨予測システム及びそれを利用した排水ポンプ運転支援システム
Pham Tracking the uncertainty in streamflow prediction through a hydrological forecasting system
Darvishi Salookolaei et al. Application of grey system theory in rainfall estimation

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: 22807312

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 22807312

Country of ref document: EP

Kind code of ref document: A1