WO2020179627A1 - パラメータ推定装置、方法、及びプログラム - Google Patents

パラメータ推定装置、方法、及びプログラム Download PDF

Info

Publication number
WO2020179627A1
WO2020179627A1 PCT/JP2020/008074 JP2020008074W WO2020179627A1 WO 2020179627 A1 WO2020179627 A1 WO 2020179627A1 JP 2020008074 W JP2020008074 W JP 2020008074W WO 2020179627 A1 WO2020179627 A1 WO 2020179627A1
Authority
WO
WIPO (PCT)
Prior art keywords
parameter
objective function
search range
input data
function value
Prior art date
Application number
PCT/JP2020/008074
Other languages
English (en)
French (fr)
Inventor
法子 横山
匡宏 幸島
達史 松林
浩之 戸田
Original Assignee
日本電信電話株式会社
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 日本電信電話株式会社 filed Critical 日本電信電話株式会社
Priority to US17/436,593 priority Critical patent/US20220171990A1/en
Publication of WO2020179627A1 publication Critical patent/WO2020179627A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
    • G06F18/2132—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on discrimination criteria, e.g. discriminant analysis
    • G06F18/21322—Rendering the within-class scatter matrix non-singular
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N5/00—Computing arrangements using knowledge-based models
    • G06N5/01—Dynamic search techniques; Heuristics; Dynamic trees; Branch-and-bound
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
    • G06F18/2135—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on approximation criteria, e.g. principal component analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N7/00—Computing arrangements based on specific mathematical models
    • G06N7/01—Probabilistic graphical models, e.g. probabilistic networks
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/21—Design or setup of recognition systems or techniques; Extraction of features in feature space; Blind source separation
    • G06F18/213—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods
    • G06F18/2132—Feature extraction, e.g. by transforming the feature space; Summarisation; Mappings, e.g. subspace methods based on discrimination criteria, e.g. discriminant analysis
    • G06F18/21322—Rendering the within-class scatter matrix non-singular
    • G06F18/21326—Rendering the within-class scatter matrix non-singular involving optimisations, e.g. using regularisation techniques

Definitions

  • the present invention relates to a parameter estimation device, a method, and a program, and more particularly to a parameter estimation device, a method, and a program for efficiently estimating parameters such as machine learning and simulation.
  • Non-Patent Document 1 Non-Patent Document 1
  • Bayesian optimization is a method of finding the optimum value (maximum value or minimum value) of the black box function based on the posterior distribution of the function obtained from the observation information.
  • the optimum value can be obtained more efficiently than a random search in which parameters are randomly selected by determining the parameters to be observed based on the probability distribution.
  • Bayesian optimization there still exists a problem that the calculation cost is high, especially when the parameters are high-dimensional or the range of possible values is large. For example, if one search takes 7 minutes, 100 times will take 10 hours or more.
  • Asynchronous parallelization of processing can be considered as a method of speeding up in Bayesian optimization.
  • Asynchronous parallel processing is a parallel processing in which each worker can start a search without waiting for the other workers to finish searching.
  • Non-Patent Document 3 there is a method for reducing the dimension of Bayesian optimization. For example, as a measure for high dimensions, there is a method of performing optimization by reducing the parameter space to low dimensions (Non-Patent Document 3). The method of Non-Patent Document 3 can be expected to find an optimal solution with a small number of observation points (the number of searches) by reducing the dimensionality.
  • Non-Patent Document 4 There is also a Bayesian optimization method for space division. For example, as a response to a high-dimensional or large-scale search section, there is a method of dividing a parameter space and optimizing each of them (Non-Patent Document 4). By narrowing the search space, it is possible to search efficiently and expect to find an optimal solution quickly.
  • the speed of parameter search can be expected.
  • the parallelization method of Non-Patent Document 2 since the past search information utilized by each worker is shared among the workers, when parallelized by M workers, a maximum of M-1 searches during the search execution are performed. There is a problem that information is not available.
  • the method of Non-Patent Document 3 only reduces the dimensionality and performs calculation, and does not perform parallelization.
  • the number of dimensions for lowering the dimension must be greater than or equal to the useful number of dimensions effective that affects the objective function.
  • the method of Non-Patent Document 4 has a problem that the process of dividing the parameter space enters every cycle, and when the time required for optimization is medium, the time required for division becomes relatively long.
  • the present invention has been made in view of the above circumstances, and an object of the present invention is to provide a parameter estimation device, method, and program capable of obtaining optimum input parameters at high speed.
  • the parameter estimation device has an input data dimension number that is the dimension number of input data, a reduction dimension number lower than the input data dimension number, and a parallel number.
  • the conversion matrix for converting the space defined by the input data dimension number into the space defined by the lower dimension number by the number of parallel number.
  • the objective function obtained by inputting the objective function value related to the observed value into a predetermined device and repeating the acquisition of the objective function value a predetermined number of times, and obtaining the objective function value from the parameter and the conversion matrix The parameter for which the objective function value is optimized based on the optimization execution unit that executes the calculation of the value a predetermined number of times in parallel and each of the objective function values obtained for each of the search ranges. And an optimum value determining unit that determines an optimum input parameter obtained from the conversion matrix.
  • the search range determination unit is based on the input data dimension number which is the dimension number of the input data, the reduction dimension number lower than the input data dimension number, and the parallel number. Then, the search for the number of parallels is performed by determining the conversion matrix for converting the space defined by the number of input data dimensions into the space defined by the number of lower dimensions for the number of parallels. The step of determining the range and the optimization executing unit are given in advance, for each of the parallel number of search ranges, the input data obtained from the parameters selected from the search range and the conversion matrix in the search range.
  • the objective function value relating to the observed value is input to a predetermined device, and the objective function value is repeated a predetermined number of times, and the objective function value is obtained from the parameter and the conversion matrix that are optimal,
  • the step of performing, in a parallel manner, repeating the determination of the function value a predetermined number of times, and the optimal value determining unit determines that the objective function value is optimal based on each of the objective function values obtained for each of the search ranges. It is characterized in that it is executed including the step of determining the optimum input parameter obtained from the above-mentioned parameters and the above-mentioned conversion matrix.
  • a program according to the present invention causes a computer to perform the parallel number based on an input data dimension number which is a dimension number of input data, a reduced dimension number lower than the input data dimension number, and a parallel number.
  • the search range for the parallel number is determined. For each of the search ranges for the number of parallels, a predetermined device that outputs the objective function value related to the observed value given in advance to the parameters selected from the search range and the input data obtained from the conversion matrix in the search range.
  • the objective function value is obtained from the parameter and the conversion matrix to be the optimum, to obtain the objective function value a predetermined number of times, Executing in parallel, based on each of the objective function values obtained for each of the search range, to determine the optimum input parameters obtained from the parameters and the transformation matrix for which the objective function value is optimal, It is a program to execute.
  • the effect that the optimum input parameters can be obtained at high speed can be obtained.
  • the asynchronous parallel processing is used to speed up the estimation.
  • the search range of the input parameter is divided by reducing the dimension of the input parameter, and the divided search range is assigned to each worker for asynchronous parallel processing.
  • the optimization may be performed in a plurality of cycles, and at that time, the search range including the optimum value calculated in the previous cycle may be preferentially assigned.
  • Bayesian optimization can be accelerated, it will be possible to make faster decisions using estimation results in the fields of tourism, transportation, medical care, etc.
  • this problem can be considered as an optimization problem that finds x in equation (1) below.
  • the conversion matrix A is a conversion matrix for converting into a space defined by the number of reduced dimensions d.
  • w is a parameter selected from the search range.
  • Non-Patent Document 3 the observed data dimensionality D, dimension reduction dimensionality d (d> useful dimensionality d effective), x ⁇ R D,
  • a ⁇ R D ⁇ d G
  • G (x ) G
  • the low-dimensional number d is not limited to a useful number of dimensions or more, and an arbitrary integer with 1 ⁇ d ⁇ D is selected as the low-dimensional number d, and the optimum solution is obtained.
  • a search can be carried out.
  • L(w) represents an objective function or an objective function value.
  • FIG. 1 is a block diagram showing an example of the configuration of the parameter estimation device 100 of this embodiment.
  • the parameter estimation device 100 having the configuration shown in FIG. 1 can be configured by a computer including a CPU, a RAM, and a ROM that stores a program for executing each processing routine described later and various data.
  • FIG. 2 is a schematic block diagram showing an example of a computer that functions as the parameter estimation device 100.
  • the parameter estimation device 100 can be realized by the computer 50 shown in FIG.
  • the computer 50 includes a CPU 51, a memory 52 as a temporary storage area, and a non-volatile storage unit 53.
  • the computer 50 also has an input / output interface (I / F) 54 to which an input / output device or the like (not shown) is connected, and a read / write (R / W) unit 55 that controls reading and writing of data to a recording medium.
  • the computer 50 includes a network I / F 56 connected to a network such as the Internet.
  • the CPU 51, the memory 52, the storage unit 53, the input / output I / F 54, the R / W unit 55, and the network I / F 56 are connected to each other via the bus 57.
  • the storage unit 53 can be realized by a Hard Disk Drive (HDD), a solid state drive (SSD), a flash memory, or the like.
  • the storage unit 53 as a storage medium stores a program for causing the computer 50 to function.
  • the CPU 51 reads the program from the storage unit 53, expands it in the memory 52, and sequentially executes the processes included in the program.
  • the parameter estimation device 100 includes a search range determination unit 102, an optimization execution unit 103, an optimum value determination unit 106, and a determination unit 107, as shown in FIG.
  • the search range determination unit 102 includes a search range determination processing unit 110 and a search range DB 111.
  • the search range determination unit 102 determines the parallel number P based on the input data dimension number D which is the dimension number of the input data, the dimension reduction dimension number d lower than the input data dimension number, and the parallel number P. Only, the search range for the parallel number P is determined by determining the conversion matrix for converting the space defined by the input data dimension number into the space defined by the lower dimension number d.
  • the optimization execution unit 103 includes area optimization execution units 104 1 to 104 N , an optimization DB 160, and an area-based optimum value DB 170.
  • Optimization execution unit 103 a search range of the parallel minutes, assigns parallel minutes space-optimized execution unit 104 1 ⁇ 104 N, by the area optimization execution unit 104 1 ⁇ 104 N in parallel the following operations Execute.
  • Each of the area optimization executing units 104 1 to 104 N inputs the input data obtained from the parameters selected from the search range and the transformation matrix into the simulator in the assigned search range, and outputs the output data and the predetermined data. Acquiring the objective function value that is the difference from the observed value is repeated a predetermined number of times.
  • each of the area optimization executing units 104 1 to 104 N approximates a function representing the relationship between the input data and the objective function value by a Gaussian process, and calculates the approximated function and the parameter for which the objective function value is optimal.
  • the acquisition function used is used to determine the parameter to be input next, input the input data obtained from the determined parameter to be input next and the conversion matrix, and obtain the objective function value a predetermined number of times. ..
  • the search range DB 111 is a database that records IDs, areas, and conversion matrices.
  • FIG. 3 is a diagram showing an example of the search range DB 111 in which data is recorded in the database table.
  • the optimization DB 160 is a database that records an ID, an area indicating a search range, a low-dimensional input parameter, an input parameter, and an objective function value.
  • FIG. 4 is a diagram showing an example of the optimization DB 160 in which data is recorded in a database table.
  • the area-specific optimum value DB170 is a database that records an ID, an area indicating a search range, an optimum low-dimensional input parameter, an optimum input parameter, and an objective function optimum value.
  • FIG. 5 is a diagram showing an example of the area-specific optimum value DB 170 in which data is recorded in a database table.
  • the optimum value determination unit 106 determines the optimum input parameters obtained from the parameters for which the objective function values are optimal and the transformation matrix, based on each of the objective function values obtained for each of the search ranges.
  • the determination unit 107 repeats each process by the search range determination unit 102, the optimization execution unit 103, and the optimum value determination unit 106 as one cycle.
  • the search range determination unit 102 determines the search range for the number of parallels by giving priority to the search range including the optimum input parameter determined in the previous cycle.
  • FIG. 6 is a flowchart showing a processing routine in parameter estimation device 100 according to the embodiment of the present invention.
  • the parameter estimation device 100 executes the processing routine shown in FIG.
  • the search range determination unit 102 performs parallel processing based on the input data dimension number D, which is the dimension number of the input data, the dimension reduction dimension number d lower than the input data dimension number, and the parallel number P.
  • the search range for the parallel number P is determined by determining the conversion matrix for converting the space defined by the input data dimension number to the space defined by the lower dimension number d by the number P. ..
  • the optimization execution unit 103 allocates the search range for the number of parallels to the region optimization execution units 104 1 to 104 N , and the region optimization execution units 104 1 to 104 N execute the next processing in parallel. To do.
  • Each of the area optimization executing units 104 1 to 104 N inputs the input data obtained from the parameters selected from the search range and the transformation matrix into the simulator in the assigned search range, and outputs the output data and the predetermined data. Acquiring the objective function value that is the difference from the observed value is repeated a predetermined number of times.
  • Each of the area optimization executing units 104 1 to 104 N approximates a function representing the relationship between the input data and the objective function value by a Gaussian process, and uses the approximated function and the parameter that optimizes the objective function value.
  • the acquisition function is used to determine the parameter to be input next, the input parameter data to be input next and the input data obtained from the transformation matrix are input to the simulator, and the objective function value is determined a predetermined number of times.
  • the optimum value determination unit 106 determines the optimum input parameter obtained from the parameter and the conversion matrix that have the optimum objective function value based on each objective function value obtained for each search range.
  • the determination unit 107 determines whether or not the predetermined condition is satisfied, and if the condition is satisfied, the process is terminated, and if the condition is not satisfied, the process returns to S100 and the process of the next cycle is executed.
  • the predetermined condition may be repeated for a predetermined number of cycles or the like. In this way, each process by the search range determination unit 102, the optimization execution unit 103, and the optimum value determination unit 106 of S100 to S104 is repeated as one cycle.
  • FIG. 7 is a flowchart showing the processing routine of the search range determining unit 102.
  • the search range determination processing unit 110 receives the input data dimension number D (observation data dimension number D), the reduced dimension number d, and the parallel number P from the user.
  • the search range determination processing unit 110 generates a designated number of parallel P search ranges.
  • the search range is, for example, an area converted into a low dimension by the method of Non-Patent Document 3. Repeat every i.
  • the transformation matrix A i ⁇ R D ⁇ d is generated as the area information for conversion to the lower dimension.
  • a method of probabilistically generating each element using an arbitrary probability distribution such as a normal distribution or a uniform distribution can be used.
  • the search range determining section 110 stores information on the generated search range (i, A i) in the search range DB 111.
  • FIG. 8 is a flowchart showing the processing routine of the optimization execution unit 103.
  • the search range generated in S100 is assigned to the optimization execution unit 103 to the area optimization execution units 104 1 to 104 N.
  • the processing routine of FIG. 8 is executed in parallel by each of the area optimization execution units 104 1 to 104 N.
  • an optimization method for example, Bayesian optimization or the like can be considered, and in the present embodiment, an example using Bayesian optimization will be described. The following describes the case where the area optimization execution unit 104 1 processes, but the same applies to the area optimization execution unit 104 N.
  • the area optimization execution unit 104 1 sets the number of random search executions to R and the number of optimization executions to B.
  • R and B the values stored in advance in the storage unit 53 or the like may be used, or the values received from the user may be used.
  • the area optimization execution unit 104 1 approximates the function f representing the relationship between the parameter w and the objective function L(w) n by a Gaussian process based on the data of the area 1 indicating the search range in the optimization DB 160. ..
  • the function f follows a Gaussian process
  • any subset follows the following Gaussian distribution.
  • K nn is the variance-covariance matrix of n ⁇ n, the variance-covariance matrix (d, d ') element k dd' kernel function k ( ⁇ , ⁇ ) using, k (w d, w It is represented by ′ d ).
  • the area optimization execution unit 104 1 determines the next observation target parameter w next that is likely to satisfy the expression (2) according to the following expression (4).
  • w next is not limited to one, and a plurality of w next may be selected in descending order of possibility.
  • ⁇ is called an acquisition function
  • several functions have been proposed (Non-Patent Document 1).
  • EI Exected Improvement
  • the acquisition function can be calculated by the following equation (5).
  • ⁇ is the average
  • ⁇ is the standard deviation
  • ⁇ and ⁇ are the CDF and PDF of the standard normal distribution, which are derived using the function f approximated by the Gaussian process in S630. Is the optimum parameter w best in the past and is acquired from the past data of the search range.
  • w next is a parameter to be input next, which is obtained by using the acquisition function. Then, (1, w, x, L(w)) is added to the optimization DB 160.
  • the area optimization execution unit 104 1 calculates w * p according to the following formula based on the data of the area 1 indicating the search range in the optimization DB 160, and calculates the area 1 of the area-specific optimum value DB 170. Update the record with (1,w * p ,x,L(w * p )).
  • the area optimization execution unit 104 1 deletes the data in the area 1 of the optimization DB 160.
  • FIG. 9 is a flowchart showing a processing routine of the optimum value determination unit 106.
  • the optimum value determination unit 106 acquires each record from the table of the area-specific optimum value DB 170.
  • the optimum value determination unit 106 calculates w * according to the following formula, and outputs a record (w * , x * , L(w)) including the optimum input parameter x * .
  • the search range determination processing unit 110 may preferentially adopt the search range including the optimum value calculated in the previous cycle. This can be dealt with by devising the method of generating the conversion matrix A or changing the definition of the objective function.
  • the parameter estimation device determines the search range for the parallel number P. For each of the search ranges corresponding to the number of parallel processes, obtaining the objective function value is repeated a predetermined number of times in parallel. As a result, Bayesian optimization can be accelerated.
  • this method can be applied to black box optimization methods other than Bayesian optimization.
  • adaptation to a genetic algorithm (GA) or evolution strategy (ES) can be considered.

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Artificial Intelligence (AREA)
  • Evolutionary Computation (AREA)
  • General Engineering & Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Mathematical Physics (AREA)
  • Software Systems (AREA)
  • Computing Systems (AREA)
  • Probability & Statistics with Applications (AREA)
  • Evolutionary Biology (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Algebra (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Computational Linguistics (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

最適な入力パラメータを高速に求めることができる。　入力データの次元数である入力データ次元数と、入力データ次元数よりも低次元の低次元化次元数と、並列数とに基づいて、並列数分だけ、入力データ次元数で定義される空間を、低次元化次元数で定義される空間に変換するための変換行列を決定することにより、並列数分の探索範囲を決定する。入力データを、シミュレータに入力し、出力データと、予め与えられた観測値との差である目的関数値を取得することを所定の回数繰り返し、獲得関数を用いて、次に入力するパラメータを決定する。決定した次に入力するパラメータ及び変換行列から得られる入力データをシミュレータに入力して、目的関数値を求めることを所定の回数繰り返すことを、並列に実行し、最適な入力パラメータを決定する。

Description

パラメータ推定装置、方法、及びプログラム
　本発明は、パラメータ推定装置、方法、及びプログラムに係り、特に、機械学習やシミュレーション等のパラメータを効率的に推定するパラメータ推定装置、方法、及びプログラムに関する。
　ＩｏＴの普及により、沢山の人やものの情報を収集できるようになってきており、取得した情報を、観光や交通、医療分野等に活かす動きが活発化している。データ解析手法として、機械学習やシミュレーション等の技術が利用されているが、これらに最適なパラメータを探索するためには、膨大な計算コストがかかることも多い。そこで、効率的なパラメータ探索の一手法として、ベイズ最適化が存在する（非特許文献１）。
　ベイズ最適化は、観測情報から得られる関数の事後分布を元に、ブラックボックス関数の最適値（最大値または最小値）を求める手法である。ベイズ最適化は、観測対象とするパラメータを、確率分布をもとに決定することで、無作為にパラメータを選定するランダムサーチなどに比べて効率的に最適値が求まることが知られている。しかし、ベイズ最適化を用いても、特にパラメータが高次元であったり、とりうる値の範囲が大きかったりする場合、計算コストがかかるという問題は依然存在している。例えば、１回の探索に７分かかる場合、１００回試行すれば１０時間以上かかることになる。
　ベイズ最適化における高速化の一手法として、処理の非同期並列化が考えられる。例えば、ベイズ最適化を並列処理する手法として、Ｍ個のワーカーにタスクを投げ、ベイズ最適化を並列で実行する研究がある（非特許文献２）。例えば１００回の探索を、Ｍ＝４個のワーカーで並列化すれば、約４倍の高速化が期待できる。この研究は、他のワーカーの待ち時間を短縮できる非同期並列処理である点も特徴である。非同期並列処理は、それぞれのワーカーが他のワーカーの探索終了を待たずに探索を開始できる並列処理である。
　また、ベイズ最適化の低次元化の手法がある。例えば、高次元への対応として、パラメータ空間を低次元に落として最適化を行う手法がある（非特許文献３）。非特許文献３の手法は、低次元化することで、少ない観測点（探索数）で最適解を見つけることが期待できる。
　また、ベイズ最適化の空間分割の手法がある。例えば、高次元や大規模な探索区間への対応として、パラメータ空間を分割し、それぞれで最適化を行う手法がある（非特許文献４）。探索空間を狭めることで、効率的に探索を行うことができ、最適解を早く見つけることが期待できる。
J.Snoek, H.Larochelle, R.P.Adams. Practical Bayesian optimization of machine learning algorithms. In Advances in Neural Information Processing Systems (NIPS), 2012. K. Kandasamy, A. Krishnamurthy, J. Schneider, B. Poczos.Parallelised Bayesian Optimisation via Thompson Sampling. In Proceedings of the Twenty-First International Conference on Articial Intelligence and Statistics (AISTATS), 2018. Z. Wang, F. Hutter, M. Zoghi, D. Matheson, N de Freitas.Bayesian Optimization in a Billion Dimensions via Random Embeddings. In Journal of Articial Intelligence Research (JAIR), 2016. Z. Wang, C. Gehring, P. Kohli, S. Jegelka. Batched large-scale bayesian optimization in high-dimensional spaces. In Proceedings of the Twenty-First International Conference on Articial Intelligence and Statistics(AISTATS), 2018.
　上記の知見を利用すれば、パラメータ探索の高速化が期待できる。しかし、非特許文献２の並列化手法では、各ワーカーが活用する過去の探索情報をワーカー間で共有するため、Ｍ個のワーカーで並列化した場合、探索実行中の最大Ｍ－１個の探索情報を利用できないという問題がある。また、非特許文献３の手法では、低次元化して計算するのみで、並列化は実施されていない。また、低次元化の次元数は、目的関数に影響を与える有用な次元数ｄｅｆｆｅｃｔｉｖｅ以上である必要がある。また、非特許文献４の手法では、パラメータ空間を分割する処理がサイクルごとに入り、最適化にかかる時間が中程度の場合、分割にかかる時間が相対的に大きくなってしまうという問題がある。
　本発明は、上記事情を鑑みて成されたものであり、最適な入力パラメータを高速に求めることができるパラメータ推定装置、方法、及びプログラムを提供することを目的とする。
　上記目的を達成するために、本発明に係るパラメータ推定装置は、入力データの次元数である入力データ次元数と、前記入力データ次元数よりも低次元の低次元化次元数と、並列数とに基づいて、前記並列数分だけ、前記入力データ次元数で定義される空間を、前記低次元化次元数で定義される空間に変換するための変換行列を決定することにより、前記並列数分の探索範囲を決定する探索範囲決定部と、前記並列数分の探索範囲の各々について、前記探索範囲において、前記探索範囲から選んだパラメータおよび前記変換行列から得られる入力データを、予め与えられた観測値に関する目的関数値を出力する所定の装置に入力し、目的関数値を取得することを所定の回数繰り返し、前記目的関数値が最適となる前記パラメータおよび前記変換行列から得られる、前記目的関数値を求めることを所定の回数繰り返すことを、並列に実行する最適化実行部と、前記探索範囲の各々について求められた目的関数値の各々に基づいて、前記目的関数値が最適となる前記パラメータ及び前記変換行列から得られる、最適な入力パラメータを決定する最適値決定部と、を含んで構成されている。
　本発明に係るパラメータ推定方法は、探索範囲決定部が、入力データの次元数である入力データ次元数と、前記入力データ次元数よりも低次元の低次元化次元数と、並列数とに基づいて、前記並列数分だけ、前記入力データ次元数で定義される空間を、前記低次元化次元数で定義される空間に変換するための変換行列を決定することにより、前記並列数分の探索範囲を決定するステップと、最適化実行部が、前記並列数分の探索範囲の各々について、前記探索範囲において、前記探索範囲から選んだパラメータおよび前記変換行列から得られる入力データを、予め与えられた観測値に関する目的関数値を出力する所定の装置に入力し、目的関数値を取得することを所定の回数繰り返し、前記目的関数値が最適となる前記パラメータおよび前記変換行列から得られる、前記目的関数値を求めることを所定の回数繰り返すことを、並列に実行するステップと、最適値決定部が、前記探索範囲の各々について求められた目的関数値の各々に基づいて、前記目的関数値が最適となる前記パラメータ及び前記変換行列から得られる、最適な入力パラメータを決定するステップと、を含んで実行することを特徴とする。
　本発明に係るプログラムは、コンピュータに、入力データの次元数である入力データ次元数と、前記入力データ次元数よりも低次元の低次元化次元数と、並列数とに基づいて、前記並列数分だけ、前記入力データ次元数で定義される空間を、前記低次元化次元数で定義される空間に変換するための変換行列を決定することにより、前記並列数分の探索範囲を決定し、前記並列数分の探索範囲の各々について、前記探索範囲において、前記探索範囲から選んだパラメータおよび前記変換行列から得られる入力データを、予め与えられた観測値に関する目的関数値を出力する所定の装置に入力し、目的関数値を取得することを所定の回数繰り返し、前記目的関数値が最適となる前記パラメータおよび前記変換行列から得られる、前記目的関数値を求めることを所定の回数繰り返すことを、並列に実行し、前記探索範囲の各々について求められた目的関数値の各々に基づいて、前記目的関数値が最適となる前記パラメータ及び前記変換行列から得られる、最適な入力パラメータを決定すること、を実行させるためのプログラムである。
　本発明のパラメータ推定装置、方法、及びプログラムによれば、最適な入力パラメータを高速に求めることができる、という効果が得られる。
本実施形態のパラメータ推定装置の構成の一例を示すブロック図である。 パラメータ推定装置として機能するコンピュータの一例を示す概略ブロック図である。 データをデータベーステーブルに記録した探索範囲ＤＢの一例を示す図である。 データをデータベーステーブルに記録した最適化ＤＢの一例を示す図である。 データをデータベーステーブルに記録した領域別最適値ＤＢの一例を示す図である。 本発明の実施の形態に係るパラメータ推定装置における処理ルーチンを示すフローチャートである。 探索範囲決定部の処理ルーチンを示すフローチャートである。 最適化実行部の処理ルーチンを示すフローチャートである。 最適値決定部の処理ルーチンを示すフローチャートである。
　以下、図面を参照して本発明の実施の形態を詳細に説明する。
　上記課題に対して、本発明の実施の形態では、複数のワーカーで処理する並列処理において、非同期並列処理を用いて推定を高速化する。
　ベイズ最適化を用いたパラメータ探索において、入力パラメータの探索範囲を、入力パラメータを低次元化することにより分割し、分割した探索範囲を各ワーカーに割り当ることで非同期並列処理を行う。
　また、最適化を複数サイクル実施しても良く、その際に、前のサイクルで算出された最適値を含む探索範囲を優先して割り当てるようにしてもよい。
　本発明の実施の形態では、非同期並列処理を用いたベイズ最適化において、利活用できない観測情報を排除するような処理をするとともに、追加計算コストがかからない方法で高速化することを目的とする。さらに、高速化手法として低次元化を用いる場合、従来手法と比較してより低次の次元数での処理を可能とすることを目的とする。
　ベイズ最適化を高速化できれば、観光や交通、医療分野等において、推定結果を用いたより早い意思決定等が可能となる。
［前提条件］
　ここで、本発明の実施の形態の前提条件について説明する。本実施形態では、局所的な場所で観察された通行人数ｙをもとに、実際の人流に近い人流シミュレーションを再現する入力パラメータｘの推定を行う問題を例に説明する。
　本問題では、入力パラメータｘおよびＧを以下のように定義する。
　ｘ：各ルートＲｊを通るエージェント数を要素とするＤ次元ベクトル
　Ｇ：ｘを入力としたシミュレーションの出力データｙｓｉｍと、実際の観測値ｙｏｂｓの差
　上記定義から、本問題は、以下（１）式のｘを求める最適化問題と考えることができる。
Figure JPOXMLDOC01-appb-M000001
                                   ・・・（１）
　さらに、探索範囲分割の一例として、入力パラメータを低次元化したものを並列数分用意して、非同期並列処理を行う。入力空間ｘに、Ｇ（ｘ）に影響を与えない探索次元や、互いに線形な関係をもつｘがある場合は、変換行列Ａを用いてＧ（ｘ）＝Ｇ（Ａｗ）となるＡとｗとが存在する。変換行列Ａは低次元化次元数ｄで定義される空間に変換するための変換行列である。ｗは探索範囲から選んだパラメータである。Ａとｗとが存在する場合は、（１）式の最適化問題を、より低次元空間での最適化問題に変換して考えることができる。非特許文献３によると、観測データ次元数Ｄ、低次元化次元数ｄ（ｄ＞有用な次元数ｄｅｆｆｅｃｔｉｖｅ）、ｘ∈ＲＤ、Ａ∈ＲＤ×ｄとすると、Ｇ（ｘ）＝Ｇ（Ａｗ）となるｗ∈Ｒｄが存在する。つまり、ｘ∈ＲＤについて、Ｇ（ｘ＊）＝Ｇ（Ａｗ＊）となるｗ＊∈Ｒｄが存在する。よって、Ｌ（ｗ）＝Ｇ（Ａｗ）とすると、本問題は、低次元空間でｗ＊を求める以下（２）式の最適化問題と考えることができる。
Figure JPOXMLDOC01-appb-M000002
                                   ・・・（２）
　ただし、本発明の実施の形態では、低次元化次元数ｄを有用な次元数以上に限定せず、１≦ｄ＜Ｄの任意の整数を低次元化次元数ｄとして選択し、最適解の探索を実施することができる。なお、Ｌ（ｗ）は、以下では、目的関数または目的関数値を表すとする。
　以上のように局所的な場所で観察された通行人数をもとに、実際の人流に近い人流シミュレーションを再現する入力パラメータの推定を行う例を用いて説明するが、機械学習やシミュレーションのパラメータ推定を行う他のデータ、事例に用いても良い。また、低次元化の手法として、非特許文献３以外の方法を用いても良い。また、並列数＝２を例に説明するが、並列数は計算環境等に応じて自由に設定して良い。さらに、探索範囲の分割手段として、入力パラメータの低次元化を複数回実行し、並列処理を行う例を説明するが、空間を複数個に分割して、それぞれの空間について並列処理を行う方法を用いても良い。
　以上の前提条件をもとに、以下、本発明の実施の形態の構成および作用について詳細に説明する。
＜本発明の実施の形態に係るパラメータ推定装置の構成＞
　次に、本発明の実施の形態に係るパラメータ推定装置の構成について説明する。
　図１は、本実施形態のパラメータ推定装置１００の構成の一例を示すブロック図である。図１に示す構成のパラメータ推定装置１００は、ＣＰＵと、ＲＡＭと、後述する各処理ルーチンを実行するためのプログラムや各種データを記憶したＲＯＭと、を含むコンピュータで構成することが出来る。
　図２は、パラメータ推定装置１００として機能するコンピュータの一例を示す概略ブロック図である。例えば、パラメータ推定装置１００は、図２に示すコンピュータ５０で実現することができる。コンピュータ５０はＣＰＵ５１、一時記憶領域としてのメモリ５２、及び不揮発性の記憶部５３を備える。また、コンピュータ５０は、入出力装置等（図示省略）が接続される入出力interface（Ｉ／Ｆ）５４、及び記録媒体に対するデータの読み込み及び書き込みを制御するread/write（Ｒ／Ｗ）部５５を備える。また、コンピュータ５０は、インターネット等のネットワークに接続されるネットワークＩ／Ｆ５６を備える。ＣＰＵ５１、メモリ５２、記憶部５３、入出力Ｉ／Ｆ５４、Ｒ／Ｗ部５５、及びネットワークＩ／Ｆ５６は、バス５７を介して互いに接続される。
　記憶部５３は、Hard Disk Drive（ＨＤＤ）、solid state drive（ＳＳＤ）、フラッシュメモリ等によって実現できる。記憶媒体としての記憶部５３には、コンピュータ５０を機能させるためのプログラムが記憶されている。ＣＰＵ５１は、プログラムを記憶部５３から読み出してメモリ５２に展開し、プログラムが有するプロセスを順次実行する。
　このパラメータ推定装置１００は、機能的には、図１に示されるように、探索範囲決定部１０２と、最適化実行部１０３と、最適値決定部１０６と、判定部１０７とを備えている。
　探索範囲決定部１０２は、探索範囲決定処理部１１０と、探索範囲ＤＢ１１１とを含む。
　探索範囲決定部１０２は、入力データの次元数である入力データ次元数Ｄと、入力データ次元数よりも低次元の低次元化次元数ｄと、並列数Ｐとに基づいて、並列数Ｐ分だけ、入力データ次元数で定義される空間を、低次元化次元数ｄで定義される空間に変換するための変換行列を決定することにより、並列数Ｐ分の探索範囲を決定する。
　最適化実行部１０３は、領域最適化実行部１０４１～１０４Ｎと、最適化ＤＢ１６０と、領域別最適値ＤＢ１７０とを含む。
　最適化実行部１０３は、並列数分の探索範囲を、並列数分の領域最適化実行部１０４１～１０４Ｎに割り当てて、領域最適化実行部１０４１～１０４Ｎによって次の処理を並列に実行する。
　領域最適化実行部１０４１～１０４Ｎの各々は、割り当てられた探索範囲において、探索範囲から選んだパラメータおよび変換行列から得られる入力データを、シミュレータに入力し、出力データと、予め与えられた観測値との差である目的関数値を取得することを所定の回数繰り返す。
　また、領域最適化実行部１０４１～１０４Ｎの各々は、入力データと目的関数値との関係を表す関数を、ガウス過程により近似し、近似した関数と、目的関数値が最適となるパラメータを用いた獲得関数を用いて、次に入力するパラメータを決定し、決定した次に入力するパラメータ及び変換行列から得られる入力データをシミュレータに入力して、目的関数値を求めることを所定の回数繰り返す。
　探索範囲ＤＢ１１１は、ＩＤ、領域、および変換行列を記録しておくデータベースである。図３は、データをデータベーステーブルに記録した探索範囲ＤＢ１１１の一例を示す図である。
　最適化ＤＢ１６０は、ＩＤ、探索範囲を示す領域、低次元入力パラメータ、入力パラメータ、および目的関数値を記録しておくデータベースである。図４は、データをデータベーステーブルに記録した最適化ＤＢ１６０の一例を示す図である。
　領域別最適値ＤＢ１７０は、ＩＤ、探索範囲を示す領域、最適低次元入力パラメータ、最適入力パラメータ、および目的関数最適値を記録しておくデータベースである。図５は、データをデータベーステーブルに記録した領域別最適値ＤＢ１７０の一例を示す図である。
　最適値決定部１０６は、探索範囲の各々について求められた目的関数値の各々に基づいて、目的関数値が最適となるパラメータ及び変換行列から得られる、最適な入力パラメータを決定する。
　判定部１０７は、探索範囲決定部１０２、最適化実行部１０３、および最適値決定部１０６による各処理を一つのサイクルとして繰り返す。サイクルを繰り返す際に、探索範囲決定部１０２は、前のサイクルで決定された最適な入力パラメータを含む探索範囲を優先して、並列数分の探索範囲を決定する。
＜本発明の実施の形態に係るパラメータ推定装置の作用＞
　次に、本発明の実施の形態に係るパラメータ推定装置１００の作用について説明する。図６は、本発明の実施の形態に係るパラメータ推定装置１００における処理ルーチンを示すフローチャートである。パラメータ推定装置１００は、図６に示す処理ルーチンを実行する。
　Ｓ１００では、探索範囲決定部１０２が、入力データの次元数である入力データ次元数Ｄと、入力データ次元数よりも低次元の低次元化次元数ｄと、並列数Ｐとに基づいて、並列数Ｐ分だけ、入力データ次元数で定義される空間を、低次元化次元数ｄで定義される空間に変換するための変換行列を決定することにより、並列数Ｐ分の探索範囲を決定する。
　Ｓ１０２では、最適化実行部１０３が、並列数分の探索範囲を、領域最適化実行部１０４１～１０４Ｎに割り当てて、領域最適化実行部１０４１～１０４Ｎによって次の処理を並列に実行する。領域最適化実行部１０４１～１０４Ｎの各々が、割り当てられた探索範囲において、探索範囲から選んだパラメータおよび変換行列から得られる入力データを、シミュレータに入力し、出力データと、予め与えられた観測値との差である目的関数値を取得することを所定の回数繰り返す。領域最適化実行部１０４１～１０４Ｎの各々が、入力データと目的関数値との関係を表す関数を、ガウス過程により近似し、近似した関数と、目的関数値が最適となるパラメータを用いた獲得関数を用いて、次に入力するパラメータを決定し、決定した次に入力するパラメータ及び変換行列から得られる入力データをシミュレータに入力して、目的関数値を求めることを所定の回数繰り返す。
　Ｓ１０４では、最適値決定部１０６が、探索範囲の各々について求められた目的関数値の各々に基づいて、目的関数値が最適となるパラメータ及び変換行列から得られる、最適な入力パラメータを決定する。
　Ｓ１０６では、判定部１０７が、所定の条件を満たすかを判定し、条件を満たす場合には処理を終了し、条件を満たさない場合にはＳ１００に戻って次のサイクルの処理を実行する。所定の条件は予め定めたサイクル数分繰り返したか等でよい。このようにして、Ｓ１００～Ｓ１０４の探索範囲決定部１０２、最適化実行部１０３、および最適値決定部１０６による各処理を一つのサイクルとして繰り返す。
　Ｓ１００の探索範囲決定部１０２の処理について説明する。図７は探索範囲決定部１０２の処理ルーチンを示すフローチャートである。
　Ｓ５１０では、探索範囲決定処理部１１０が、ユーザから、入力データ次元数Ｄ（観測データ次元数Ｄ）、低次元化次元数ｄ、および並列数Ｐを受け付ける。
　Ｓ５２０の処理ステップは、ｉ＝０から実行ごとにｉ＝ｉ＋１でカウントアップし、ｉ＜Ｐの条件で繰り返し、並列数Ｐ回分実行する。
　Ｓ５２０では、探索範囲決定処理部１１０は、指定された並列数Ｐ個の探索範囲を生成する。探索範囲は、例えば非特許文献３の方法で低次元に変換した領域とする。ｉごとに繰り返す。具体的には、低次元に変換するための領域情報として、変換行列Ａｉ∈ＲＤ×ｄを生成する。変換行列Ａの生成には、各要素を正規分布や一様分布など任意の確率分布を用いて確率的に生成する方法などが利用できる。
　Ｓ５３０では、探索範囲決定処理部１１０が、生成した探索範囲に関する情報（ｉ，Ａｉ）を探索範囲ＤＢ１１１に格納する。
　次にＳ１０２の最適化実行部１０３の処理について説明する。図８は、最適化実行部１０３の処理ルーチンを示すフローチャートである。最適化実行部１０３には、Ｓ１００で生成された探索範囲が、領域最適化実行部１０４１～１０４Ｎに割り当てられる。本実施形態ではＮ＝２である場合を例に説明する。図８の処理ルーチンは、領域最適化実行部１０４１～１０４Ｎの各々でそれぞれ並列に実行する。また、最適化手法としては、例えばベイズ最適化等を利用することが考えられ、本実施形態ではベイズ最適化を用いた例を説明する。以下は、領域最適化実行部１０４１が処理する場合を例に説明するが、領域最適化実行部１０４Ｎでも同様である。
　Ｓ６１０では、領域最適化実行部１０４１が、ランダムサーチ実行回数をＲ、最適化実行回数をＢとする。ＲおよびＢの値は、あらかじめ記憶部５３等に保存した値を用いても良いし、ユーザから入力を受け付けた値を用いても良い。
　Ｓ６２０の処理ステップは、ｉ＝０から実行ごとにｉ＝ｉ＋１でカウントアップし、ｉ＜Ｒの条件で繰り返し、ランダムサーチ実行回数のＲ回分実行する。
　Ｓ６２０では、領域最適化実行部１０４１が、当該探索範囲を示す領域１の中でランダムにパラメータｗを選び、ｘ＝Ａ１ｗとしてシミュレータに入力し、目的関数値Ｌ（ｗ）を取得する。そして、最適化ＤＢ１６０に（１，ｗ，ｘ，Ｌ（ｗ））を追加する。
　以下のＳ６３０～Ｓ６５０の処理ステップは、ｉ＝０から実行ごとにｉ＝ｉ＋１でカウントアップし、ｉ＜Ｂの条件で繰り返し、最適化実行回数のＢ回分実行する。
　Ｓ６３０では、領域最適化実行部１０４１が、最適化ＤＢ１６０の当該探索範囲を示す領域１のデータに基づき、パラメータｗと目的関数Ｌ（ｗ）ｎの関係を表す関数ｆをガウス過程により近似する。ここで、関数ｆがガウス過程に従うとすると、任意の部分集合は次のガウス分布に従う。
Ｐ（ｆ）＝Ｎ（ｆ｜０，Ｋｎｎ）
　ただし、Ｋｎｎはｎ×ｎの分散共分散行列であり、分散共分散行列の（ｄ，ｄ′）要素ｋｄｄ′はカーネル関数ｋ（・，・）を用いて、ｋ（ｗｄ，ｗ′ｄ）で表される。
　そして、ｙ＝Ｌ（ｗ）が与えられたもとでの未知の入力変数ｗ＊の出力変数ｙ＊の事後確率は、以下（３）式のガウス分布で与えられる。
Figure JPOXMLDOC01-appb-M000003
                                                 ・・・（３）
　ただし、Ｃｎｎ＝Ｋｎｎ＋σ２Ｉｎと定義し、Ｉｎはｎ×ｎの単位行列を表す。よって、（３）式を用いて予測ができる。
　Ｓ６４０では、領域最適化実行部１０４１が、以下（４）式に従い、（２）式を満たす可能性が高い次の観測対象パラメータｗｎｅｘｔを決定する。この時、ｗｎｅｘｔは一つではなくて、可能性が高い順に複数選んでもよい。
Figure JPOXMLDOC01-appb-M000004
                            ・・・（４）
　ここで、αは獲得関数と呼ばれ、いくつかの関数が提案されている（非特許文献１）。例えば、ＥＩ（Ｅｘｐｅｃｔｅｄ　Ｉｍｐｒｏｖｅｍｅｎｔ）を用いる場合、獲得関数は以下（５）式で計算できる。
Figure JPOXMLDOC01-appb-M000005
                                                 ・・・（５）
　ここで、μは平均、σは標準偏差、Φおよびφは標準正規分布のＣＤＦおよびＰＤＦであり、Ｓ６３０でガウス過程により近似した関数ｆを用いて導出される。
Figure JPOXMLDOC01-appb-I000006

は過去の最適なパラメータｗｂｅｓｔであり、当該探索範囲の過去のデータから取得する。
　Ｓ６５０では、領域最適化実行部１０４１が、ｗ＝ｗｎｅｘｔにＡをかけて求めた入力パラメータｘをシミュレータに入力し、目的関数値Ｌ（ｗ）を取得する。ｗｎｅｘｔは、獲得関数を用いて求めた、次に入力するパラメータである。そして、最適化ＤＢ１６０に（１，ｗ，ｘ，Ｌ（ｗ））を追加する。
　Ｓ６６０では、領域最適化実行部１０４１が、最適化ＤＢ１６０の当該探索範囲を示す領域１のデータをもとに、以下の式に従ってｗ＊ ｐを算出し、領域別最適値ＤＢ１７０の領域１のレコードを（１，ｗ＊ ｐ，ｘ，Ｌ（ｗ＊ ｐ））で更新する。
Figure JPOXMLDOC01-appb-I000007

　Ｓ６７０では、領域最適化実行部１０４１が、最適化ＤＢ１６０の領域１のデータを削除する。
　次にＳ１０４の最適値決定部１０６の処理について説明する。図９は、最適値決定部１０６の処理ルーチンを示すフローチャートである。
　Ｓ７１０では、最適値決定部１０６が、領域別最適値ＤＢ１７０のテーブルから各レコードを取得する。
　Ｓ７２０では、最適値決定部１０６が、以下の式に従ってｗ＊を算出し、最適な入力パラメータｘ＊を含むレコード（ｗ＊，ｘ＊，Ｌ（ｗ））を出力する。
Figure JPOXMLDOC01-appb-I000008

　判定部１０７により複数サイクルで繰り返す場合、探索範囲決定処理部１１０で、前のサイクルで算出された最適値を含む探索範囲を優先して採用しても良い。これは、変換行列Ａの生成の方法を工夫する、または目的関数の定義を変えることで対応できる。目的関数の定義を変える場合、例えばある点ｘｐｒｉｏｒの周辺に最適点が存在することが分かるならば、目的関数をＬ（ｗ）＝Ｇ（Ａｗ＋ｘｐｒｉｏｒ）と定義する。このように目的関数を定義することでＬ（ｗ）により求まる点を通る平面での最適化を考えることができる。このように、判定部１０７の繰り返しにおいて、前のサイクルで決定された最適な入力パラメータｘｐｒｉｏｒを含む探索範囲を優先して採用するように、次のサイクルの目的関数を定義するようにしても良い。
　以上説明したように、本発明の実施の形態に係るパラメータ推定装置によれば、並列数Ｐ分の探索範囲を決定する。並列数分の探索範囲の各々について、目的関数値を求めることを所定の回数繰り返すことを、並列に実行する。これにより、ベイズ最適化を高速化できる。
　また、探索範囲を分割し小さくすることで、ベイズ最適化の高速化を実現する。また、探索範囲を低次元化により分割する処理とすることで、利活用できない観測情報を排除した処理が実現できる。また、簡易な方法で探索範囲を分割することにより、最適化にかかる時間が中程度の場合でも、分割にかかる時間が相対的に大きくならないようにすることができる。また、探索範囲を分割することで、ベイズ最適化が適応困難な高次元や探索範囲の広いパラメータ探索を可能とする。また、探索範囲の分割として低次元化を用いる場合、低次元化された空間の探索を並列処理、および複数サイクル実施することにより、１≦ｄ＜Ｄの任意の整数を低次元化次元数ｄとして選択することができ、より低次元での最適解の探索を可能とする。
　なお、本発明は、上述した実施の形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
　たとえば、本手法は、ベイズ最適化以外のブラックボックス最適化手法にも適用可能である。例えば、遺伝的アルゴリズム（ＧＡ）や進化戦略（ＥＳ）などへの適応が考えられる。
　上述した実施形態では、ベイズ最適化において、ガウス過程による近似を用いる場合を例に説明したが、ランダムフォレストなど、別の確率モデルを用いた近似手法を用いてもよい。
　上述した実施形態では、ベイズ最適化に適用した場合に、所定の装置をシミュレータとして、入力パラメータをシミュレータに入力して目的関数値を得る場合を例に説明したが、各々の最適化手法に対応した所定の装置に入力し、目的関数値を取得すればよい。
１００ パラメータ推定装置
１０２ 探索範囲決定部
１０３ 最適化実行部
１０４１～Ｎ 領域最適化実行部
１０６ 最適値決定部
１０７ 判定部
１１０ 探索範囲決定処理部

Claims (7)

  1. 　入力データの次元数である入力データ次元数と、前記入力データ次元数よりも低次元の低次元化次元数と、並列数とに基づいて、前記並列数分だけ、前記入力データ次元数で定義される空間を、前記低次元化次元数で定義される空間に変換するための変換行列を決定することにより、前記並列数分の探索範囲を決定する探索範囲決定部と、
    　前記並列数分の探索範囲の各々について、
    　前記探索範囲において、前記探索範囲から選んだパラメータおよび前記変換行列から得られる入力データを、予め与えられた観測値に関する目的関数値を出力する所定の装置に入力し、目的関数値を取得することを所定の回数繰り返し、
    　前記目的関数値が最適となる前記パラメータおよび前記変換行列から得られる、前記目的関数値を求めることを所定の回数繰り返すことを、並列に実行する最適化実行部と、
    　前記探索範囲の各々について求められた目的関数値の各々に基づいて、前記目的関数値が最適となる前記パラメータ及び前記変換行列から得られる、最適な入力パラメータを決定する最適値決定部と、
    　を含むパラメータ推定装置。
  2. 　前記最適化実行部は、前記並列数分の探索範囲の各々について、前記目的関数値を取得したあとに、
    　前記入力データと前記目的関数値との関係を表す関数を、確率モデルにより近似し、前記近似した関数と、前記目的関数値が最適となる前記パラメータを用いた獲得関数を用いて、次に入力するパラメータを決定し、決定した前記次に入力するパラメータおよび前記変換行列から得られる入力データを、前記所定の装置に入力し、前記目的関数値を求めることを所定の回数繰り返すことを、並列に実行する請求項１に記載のパラメータ推定装置。
  3. 　前記最適化実行部は、前記並列数分の探索範囲の各々について、
    　前記探索範囲において、前記探索範囲から選んだパラメータおよび前記変換行列から得られる入力データを、シミュレータに入力し、出力データと、前記目的関数値を取得することを所定の回数繰り返し、
    　前記獲得関数を用いて、次に入力するパラメータを決定し、決定した前記次に入力するパラメータおよび前記変換行列から得られる入力データを、前記シミュレータに入力し、前記目的関数値を求めることを所定の回数繰り返すことを、並列に実行する請求項２に記載のパラメータ推定装置。
  4. 　前記探索範囲決定部、前記最適化実行部、および前記最適値決定部による各処理を一つのサイクルとして繰り返す判定部を更に含む請求項１～請求項３の何れか１項に記載のパラメータ推定装置。
  5. 　前記繰り返しにおいて、前記探索範囲決定部は、前のサイクルで決定された前記最適な入力パラメータを含む探索範囲を優先して、前記並列数分の探索範囲を決定する請求項４に記載のパラメータ推定装置。
  6. 　探索範囲決定部が、入力データの次元数である入力データ次元数と、前記入力データ次元数よりも低次元の低次元化次元数と、並列数とに基づいて、前記並列数分だけ、前記入力データ次元数で定義される空間を、前記低次元化次元数で定義される空間に変換するための変換行列を決定することにより、前記並列数分の探索範囲を決定するステップと、
    　最適化実行部が、前記並列数分の探索範囲の各々について、
    　前記探索範囲において、前記探索範囲から選んだパラメータおよび前記変換行列から得られる入力データを、予め与えられた観測値に関する目的関数値を出力する所定の装置に入力し、目的関数値を取得することを所定の回数繰り返し、
    　前記目的関数値が最適となる前記パラメータおよび前記変換行列から得られる、前記目的関数値を求めることを所定の回数繰り返すことを、並列に実行するステップと、
    　最適値決定部が、前記探索範囲の各々について求められた目的関数値の各々に基づいて、前記目的関数値が最適となる前記パラメータ及び前記変換行列から得られる、最適な入力パラメータを決定するステップと、
    　を含むパラメータ推定方法。
  7. 　コンピュータに、
    　入力データの次元数である入力データ次元数と、前記入力データ次元数よりも低次元の低次元化次元数と、並列数とに基づいて、前記並列数分だけ、前記入力データ次元数で定義される空間を、前記低次元化次元数で定義される空間に変換するための変換行列を決定することにより、前記並列数分の探索範囲を決定し、
    　前記並列数分の探索範囲の各々について、
    　前記探索範囲において、前記探索範囲から選んだパラメータおよび前記変換行列から得られる入力データを、予め与えられた観測値に関する目的関数値を出力する所定の装置に入力し、目的関数値を取得することを所定の回数繰り返し、
    　前記目的関数値が最適となる前記パラメータおよび前記変換行列から得られる、前記目的関数値を求めることを所定の回数繰り返すことを、並列に実行し、
    　前記探索範囲の各々について求められた目的関数値の各々に基づいて、前記目的関数値が最適となる前記パラメータ及び前記変換行列から得られる、最適な入力パラメータを決定すること、
    　を実行させるためのプログラム。
PCT/JP2020/008074 2019-03-05 2020-02-27 パラメータ推定装置、方法、及びプログラム WO2020179627A1 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US17/436,593 US20220171990A1 (en) 2019-03-05 2020-02-27 Parameter estimation device, method and program

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2019-039694 2019-03-05
JP2019039694A JP7124768B2 (ja) 2019-03-05 2019-03-05 パラメータ推定装置、方法、及びプログラム

Publications (1)

Publication Number Publication Date
WO2020179627A1 true WO2020179627A1 (ja) 2020-09-10

Family

ID=72337740

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2020/008074 WO2020179627A1 (ja) 2019-03-05 2020-02-27 パラメータ推定装置、方法、及びプログラム

Country Status (3)

Country Link
US (1) US20220171990A1 (ja)
JP (1) JP7124768B2 (ja)
WO (1) WO2020179627A1 (ja)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2022086962A (ja) * 2020-11-30 2022-06-09 財團法人工業技術研究院 適応的探索空間を用いた、パラメータ最適化のためのシステムおよび方法、ならびにそれを用いたユーザインターフェース
US11960253B2 (en) 2020-11-30 2024-04-16 Industrial Technology Research Institute System and method for parameter optimization with adaptive search space and user interface using the same

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2016523402A (ja) * 2013-05-30 2016-08-08 プレジデント アンド フェローズ オブ ハーバード カレッジ ベイズの最適化を実施するためのシステムおよび方法
CN109376869A (zh) * 2018-12-25 2019-02-22 中国科学院软件研究所 一种基于异步贝叶斯优化的机器学习超参优化系统及方法

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2016523402A (ja) * 2013-05-30 2016-08-08 プレジデント アンド フェローズ オブ ハーバード カレッジ ベイズの最適化を実施するためのシステムおよび方法
CN109376869A (zh) * 2018-12-25 2019-02-22 中国科学院软件研究所 一种基于异步贝叶斯优化的机器学习超参优化系统及方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
QIAN, HONG ET AL.: "Scaling Simultaneous Optimistic Optimization for High- Dimensional Non-Convex Functions with Low Effective Dimensions", PROCEEDINGS OF THE THIRTIETH AAAI CONFERENCE ON ARTIFICIAL INTELLIGENCE (AAAI-16), 12 February 2016 (2016-02-12) - 17 February 2016 (2016-02-17), pages 2000 - 2006, XP055736679 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2022086962A (ja) * 2020-11-30 2022-06-09 財團法人工業技術研究院 適応的探索空間を用いた、パラメータ最適化のためのシステムおよび方法、ならびにそれを用いたユーザインターフェース
JP7158512B2 (ja) 2020-11-30 2022-10-21 財團法人工業技術研究院 適応的探索空間を用いた、パラメータ最適化のためのシステムおよび方法、ならびにそれを用いたユーザインターフェース
US11960253B2 (en) 2020-11-30 2024-04-16 Industrial Technology Research Institute System and method for parameter optimization with adaptive search space and user interface using the same

Also Published As

Publication number Publication date
JP7124768B2 (ja) 2022-08-24
JP2020144530A (ja) 2020-09-10
US20220171990A1 (en) 2022-06-02

Similar Documents

Publication Publication Date Title
Lei et al. GCN-GAN: A non-linear temporal link prediction model for weighted dynamic networks
Hong et al. Dropnas: Grouped operation dropout for differentiable architecture search
Chen et al. Bridging the gap between stochastic gradient MCMC and stochastic optimization
WO2020179627A1 (ja) パラメータ推定装置、方法、及びプログラム
WO2010064414A1 (ja) 遺伝子クラスタリングプログラム、遺伝子クラスタリング方法及び遺伝子クラスター解析装置
US20060217939A1 (en) Time series analysis system, time series analysis method, and time series analysis program
CN111258767A (zh) 复杂系统仿真应用的云计算资源智能分配方法与装置
Krishnakumar et al. Nas-bench-suite-zero: Accelerating research on zero cost proxies
Xu et al. Continuous-action reinforcement learning with fast policy search and adaptive basis function selection
Lin et al. Aacp: Model compression by accurate and automatic channel pruning
JP2022070386A (ja) 学習方法、系列解析方法、学習装置、系列解析装置、及びプログラム
JP2018067227A (ja) データ分析装置、データ分析方法、データ分析処理プログラム
Garant et al. Accelerating multi-agent reinforcement learning with dynamic co-learning
JP2021524112A (ja) 情報処理装置、制御方法および非一時的記憶媒体
JP6398991B2 (ja) モデル推定装置、方法およびプログラム
JP2005135010A (ja) プロセスの状態類似事例検索方法および状態予測方法並びに記憶媒体
WO2020218246A1 (ja) 最適化装置、最適化方法、及びプログラム
Tomczak et al. Policy optimization through approximate importance sampling
Harland et al. Path ensembles and path sampling in nonequilibrium stochastic systems
JP7310884B2 (ja) パラメータ推定装置、パラメータ推定方法、及びパラメータ推定プログラム
Du et al. Novel quantization strategies for linear prediction with guarantees
Gruner et al. Pseudo-Likelihood Inference
JPH06149866A (ja) 解探索装置
Cao et al. Efficient model performance estimation via feature histories
CN116993181A (zh) 基于rbf-arx的综合能源系统时序概率多元负荷预测方法及系统

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

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

Country of ref document: EP

Kind code of ref document: A1