WO2021220399A1 - 推定装置、推定システム、推定方法及びプログラム - Google Patents

推定装置、推定システム、推定方法及びプログラム Download PDF

Info

Publication number
WO2021220399A1
WO2021220399A1 PCT/JP2020/018115 JP2020018115W WO2021220399A1 WO 2021220399 A1 WO2021220399 A1 WO 2021220399A1 JP 2020018115 W JP2020018115 W JP 2020018115W WO 2021220399 A1 WO2021220399 A1 WO 2021220399A1
Authority
WO
WIPO (PCT)
Prior art keywords
function
estimation
rate
kernel
spatial region
Prior art date
Application number
PCT/JP2020/018115
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/919,910 priority Critical patent/US20230205841A1/en
Priority to JP2022518485A priority patent/JP7392838B2/ja
Priority to PCT/JP2020/018115 priority patent/WO2021220399A1/ja
Publication of WO2021220399A1 publication Critical patent/WO2021220399A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/17—Function evaluation by approximation methods, e.g. inter- or extrapolation, smoothing, least mean square method
    • G06F17/175—Function evaluation by approximation methods, e.g. inter- or extrapolation, smoothing, least mean square method of multidimensional data
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/18—Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06N—COMPUTING ARRANGEMENTS BASED ON SPECIFIC COMPUTATIONAL MODELS
    • G06N99/00—Subject matter not provided for in other groups of this subclass

Definitions

  • the present invention relates to an estimation device, an estimation system, an estimation method and a program.
  • Gaussian Cox Process is known as a technique for estimating the event occurrence rate (hereinafter referred to as "rate function") at each point in the space using data representing events that occur in one-dimensional or multidimensional space.
  • rate function the event occurrence rate
  • Non-Patent Document 1 proposes an approximate method for estimating a rate function using the variational Bayesian method.
  • the approximation using the variational Bayesian method may cause errors or biases in the estimation results, and there is a possibility that an estimated value that is significantly different from the true rate function will be output.
  • it is difficult to quantitatively evaluate the magnitude of the error of the estimated value and there is a possibility that an incorrect estimated value is adopted as the correct one.
  • One embodiment of the present invention has been made in view of the above points, and an object thereof is to analytically estimate a rate function.
  • the estimation device includes a set of observation points representing events observed in a spatial region of a predetermined dimension, an number of observations representing the number of times the observation is made, and a predetermined number of observations. It has a rate function estimation unit that analytically estimates a rate function for obtaining the occurrence rate of the event in the space region by inputting a first function satisfying the condition and parameters of the first function. It is characterized by that.
  • the rate function can be estimated analytically.
  • FIG. It is a figure which shows an example of the functional structure of the estimation apparatus in Example 1.
  • FIG. It is a flowchart which shows an example of the flow of the estimation process in Example 1.
  • FIG. It is a figure which shows an example of the functional structure of the estimation apparatus in Example 2.
  • FIG. It is a flowchart which shows an example of the flow of the estimation process in Example 2.
  • FIG. It is a flowchart which shows an example of the flow of the estimation process in Example 3.
  • FIG. It is a figure which shows an example of the hardware composition of the estimation apparatus which concerns on this embodiment.
  • the estimation device 10 capable of analytically estimating the rate function based on the Gaussian Cox Process (hereinafter, abbreviated as “GCP”) will be described.
  • GCP Gaussian Cox Process
  • FIG. 1 is a diagram showing an example of the functional configuration of the estimation device 10 in the first embodiment.
  • the estimation device 10 in the first embodiment has an estimator generation unit 101 and an estimation unit 102. Further, it is assumed that the estimation device 10 in the first embodiment is given an estimation target point 1000, point event data 1100, a GCP kernel 1200, and a GCP kernel parameter 1300.
  • the estimation target point 1000, the point event data 1100, the GCP kernel 1200, and the GCP kernel parameter 1300 may be stored in a storage device such as an HDD (Hard Disk Drive) or an SSD (Solid State Drive), for example. It may be input by a user's operation or the like, or may be received from another device connected via a communication network.
  • HDD Hard Disk Drive
  • SSD Solid State Drive
  • the estimation target point 1000 is a set ⁇ s 1 , s 2 , ..., S P ⁇ of P points (coordinate points) for which the value of the rate function is to be estimated.
  • the point event data 1100 is a set of the spatial region T where the observation was performed, the number of observations R, and the points (coordinate points) of N events observed in the spatial region T ⁇ t 1 , t 2 , ... ⁇ ⁇ , T N ⁇ .
  • the spatial region T is a predetermined spatial region of an arbitrary dimension (one-dimensional or multidimensional).
  • the GCP kernel 1200 is a function h (t, t') that satisfies the following conditions 1 and 2. However, t and t'are arbitrary points in the space region T, respectively.
  • ⁇ (t) is an eigenfunction
  • GCP kernel the function satisfying the above conditions 1 and 2 will be referred to as a "GCP kernel” or a "GCP kernel function”.
  • GCP kernel function An example of a GCP kernel is shown below.
  • Example of GCP kernel 1 The function h shown in the following equation (2) defined by the semi-fixed value kernel k (t, t') and the upper limit of its eigenvalue k max is the GCP kernel.
  • k max is evaluated as the maximum value of the function obtained by integrating the semi-fixed value kernel k (t, t') with respect to one variable (t'in the equation (3)) as shown in the following equation (3). Will be done.
  • the semi-normal kernel k (t, t') includes the Gaussian kernel k Gauss (t, t') having the intensity ⁇ and the covariance matrix ⁇ as parameters, and the Matern kernel k Matern having the intensity ⁇ and the case ⁇ as parameters. There are (t, t') and the like.
  • the Gaussian kernel k Gauss (t, t') is shown in the following equation (4)
  • the Matern kernel k Matern (t, t') is shown in the following equation (5).
  • ⁇ represents a transpose operation on vectors and matrices.
  • the Wiener kernel k Wiener (t, t') can also be used.
  • Example of GCP kernel 2 The solution of the second kind Fredholm type integral equation defined using the semi-fixed value kernel k (t, t') is the GCP kernel. That is, the solution h (t, t') of the integral equation shown in the following equation (6) is the GCP kernel.
  • the GCP kernel parameter 1300 is a parameter value of the GCP kernel h (t, t').
  • the number of parameters of the GCP kernel h (t, t') and the range of values that the parameters can take depend on the definition of the GCP kernel h (t, t').
  • the Gaussian kernel k Gauss (t, t') is used for the semi-definite matrix k (t, t') and the function h shown in the above equation (2) is the GCP kernel
  • the parameter is the intensity ⁇ .
  • the covariance matrix ⁇ is the covariance matrix ⁇ .
  • the GCP kernel parameter 1300 includes a parameter ⁇ representing the square root of the average value of the rate function.
  • the parameter ⁇ shall take a non-negative value.
  • each parameter included in the GCP kernel parameter 1300 is usually determined based on various prior knowledge, restrictions, and the like for each case in which the estimation device 10 according to the present embodiment is applied. For example, for an event in which a user of an EC site purchases a product with a strong seasonality, the Matern kernel k Matern (t, t') is used to convert the function h shown in the above equation (2) to GCP. In the case of a kernel, the parameter ⁇ can be considered to be a value representing several weeks. Further, for example, it is conceivable to estimate each rate function using some candidate parameter values, compare the estimation results, and adopt a parameter value that is subjectively considered to be valid.
  • the estimator generator 101 inputs the point event data 1100, the GCP kernel 1200, and the GCP kernel parameter 1300 given to the estimator 10, and generates and outputs the rate function estimator 1400.
  • the rate function estimator 1400 is a function ⁇ (t) for obtaining an estimated value of the rate function at an arbitrary point t in the space region T.
  • This function ⁇ (t) is expressed by the following equation (7).
  • a vector composed of values, z is a vector composed of the reciprocals of the square roots of the estimated values of the rate functions on the coordinate points t 1 , t 2 , ..., T N of N events. That is,
  • the rate function estimator 1400 (that is, ⁇ (t)) is
  • the estimator generator 101 outputs the rate function estimator 1400 by outputting these five elements.
  • the estimator generation unit 101 generates (calculates) z, which is an unknown number among these five elements, by calculation, and then outputs the rate function estimator 1400.
  • This conditional equation is an N-element simultaneous quadratic equation shown in the following equation (8).
  • n' is also 1 ⁇ n' ⁇ N.
  • the estimation unit 102 inputs the estimation target point 1000 given to the estimation device 10 and the rate function estimator 1400 output by the estimator generation unit 101, and calculates and outputs the estimated value 1500. That is, the estimation unit 102 substitutes each of the P points s 1 , s 2 , ..., And s P into ⁇ (t) shown in the above equation (7) to obtain the estimated value of the rate function.
  • the set ⁇ (s 1 ), ⁇ (s 2 ), ..., ⁇ (s P ) ⁇ is calculated. This ⁇ (s 1 ), ⁇ (s 2 ), ..., ⁇ (s P ) ⁇ is an estimated value of 1500.
  • FIG. 2 is a flowchart showing an example of the flow of the estimation process in the first embodiment.
  • Step S101 First, the estimator generator 101 inputs the given point event data 1100, the GCP kernel 1200, and the GCP kernel parameter 1300, and generates and outputs the rate function estimator 1400.
  • Step S102 Then, the estimation unit 102 inputs the given estimation target point 1000 and the rate function estimator 1400 output in step S101 above, and calculates and outputs the estimated value 1500. As a result, a set of estimated values of the rate function ⁇ (s 1 ), ⁇ (s 2 ), ..., ⁇ (s P ) ⁇ is obtained. At this time, since the estimation device 10 in the first embodiment estimates the rate function analytically based on GCP without using approximation, it is possible to prevent a situation in which an unintended error or bias is included in the estimation result. Can be done.
  • the process of step S101 may be executed in advance before the process of step S102 is executed (in this case, the rate function estimator 1400 is stored in a storage device or the like, and the process of step S101 is performed. In the process, the rate function estimator 1400 may be read from a storage device or the like.) At this time, when the process of step S102 is executed a plurality of times using the same rate function estimator 1400, the process of step S101 may be executed once. Further, for example, when the process of the above step S101 is executed by another device and the rate function estimator 1400 is given from the other device, the estimator 10 executes only the process of the above step S102. (In this case, the system is composed of the other device and the estimator 10 and the other device has the estimator generator 101, while the estimator 10 does not have the estimator generator 101. May be.).
  • Example 2 Next, Example 2 of this embodiment will be described.
  • the value of each parameter included in the GCP kernel parameter 1300 is not determined subjectively but is determined without subjectivity will be described.
  • FIG. 3 is a diagram showing an example of the functional configuration of the estimation device 10 in the second embodiment.
  • the estimation device 10 in the second embodiment has an estimator generation unit 101, an estimation unit 102, and a parameter determination unit 103. Further, it is assumed that the estimation device 10 in the second embodiment is provided with the estimation target point 1000, the point event data 1100, and the GCP kernel 1200.
  • the estimated target points 1000, the point event data 1100, and the GCP kernel 1200 may be stored in a storage device such as an HDD or SSD, may be input by a user's operation, or may be input to a communication network. It may be received from other devices connected via.
  • the estimator generation unit 101 inputs the GCP kernel parameter 1300 output by the parameter determination unit 103.
  • the parameter determination unit 103 inputs the point event data 1100 and the GCP kernel 1200 given to the estimation device 10, and outputs the GCP kernel parameter 1300. That is, the parameter determination unit 103 sets the spatial region T, the number of observations R, the set of N event points ⁇ t 1 , t 2 , ..., T N ⁇ and the GCP kernel h (t, t'). By cross-validation, the value of each parameter included in the GCP kernel parameter 1300 is determined to be the optimum value.
  • the parameter determination unit 103 divides, for example, a set of points of N events ⁇ t 1 , t 2 , ..., T N ⁇ into K batches (subsets), and K-
  • K may be determined according to the spatial region T and the number of observations R.
  • FIG. 4 is a flowchart showing an example of the flow of the estimation process in the second embodiment.
  • Step S201 First, the parameter determination unit 103 inputs the given point event data 1100 and the GCP kernel 1200, and outputs the GCP kernel parameter 1300.
  • Step S202 Next, the estimator generator 101 inputs the given point event data 1100 and the GCP kernel 1200 and the GCP kernel parameter 1300 output in step S201 above to generate the rate function estimator 1400. Generate and output.
  • Step S203 Then, the estimation unit 102 inputs the given estimation target point 1000 and the rate function estimator 1400 output in step S202 above, and calculates and outputs the estimated value 1500. As a result, a set of estimated values of the rate function ⁇ (s 1 ), ⁇ (s 2 ), ..., ⁇ (s P ) ⁇ is obtained. Moreover, since the estimation device 10 in the second embodiment can optimize the GCP kernel parameter 1300, for example, the knowledge and experience regarding the prior knowledge and restrictions of the case where the estimation device 10 according to the present embodiment is applied are insufficient. Even a user can obtain an accurate estimate of the rate function.
  • steps S201 and S202 may be executed in advance, or the processes of steps S201 and S202 may be executed by another device.
  • Example 3 Next, Example 3 of this embodiment will be described. In the third embodiment, a case where the error of the estimated value of the rate function is also calculated will be described.
  • FIG. 5 is a diagram showing an example of the functional configuration of the estimation device 10 in the third embodiment.
  • the estimation device 10 in the third embodiment has an estimator generation unit 101, an estimation unit 102, and an error calculation unit 104. Further, it is assumed that the estimation device 10 in the third embodiment is given an estimation target point 1000, point event data 1100, a GCP kernel 1200, and a GCP kernel parameter 1300.
  • the estimation target point 1000, the point event data 1100, the GCP kernel 1200, and the GCP kernel parameter 1300 may be stored in a storage device such as an HDD or SSD, or may be input by a user operation or the like. However, it may be received from another device connected via a communication network.
  • estimator generator 101 and the estimator 102 are the same as those in the first embodiment, the description thereof will be omitted.
  • the error calculation unit 104 inputs the estimation target point 1000 given to the estimation device 10, the rate function estimator 1400 output by the estimator generation unit 101, and the estimation value 1500 output by the estimation unit 102. , Calculate and output the error 1600.
  • the error 1600 is an estimation error (covariance matrix) at P points s 1 , s 2 , ..., S P to be estimated with respect to the square root of the rate function.
  • These estimation errors ⁇ (t) are calculated by the following equation (9).
  • FIG. 6 is a flowchart showing an example of the flow of the estimation process in the third embodiment.
  • Step S301 First, the estimator generator 101 inputs the given point event data 1100, the GCP kernel 1200, and the GCP kernel parameter 1300, and generates and outputs the rate function estimator 1400.
  • Step S302 Next, the estimation unit 102 inputs the given estimation target point 1000 and the rate function estimator 1400 output in step S301 above, and calculates and outputs the estimated value 1500.
  • Step S303 Then, the error calculation unit 104 inputs the given estimation target point 1000, the rate function estimator 1400 output in the above step S301, and the estimated value 1500 output in the above step S302. Therefore, the error 1600 is calculated and output.
  • the estimation device 10 in the third embodiment can quantitatively evaluate the magnitude of the error of the estimated value of the rate function.
  • step S301 may be executed in advance, or the process of step S301 may be executed by another device.
  • the rate function by giving the point event data 1100 including the purchase time data of one user, the spatial area T, and the number of observations R, and the GCP kernel 1200 and the GCP kernel parameter 1300 selected by the person in charge of advertising strategy to the estimation device 10.
  • the estimator 1400 is output.
  • the person in charge of the advertising strategy sets the time point sequence ⁇ s 1 , s 2 , ..., S P ⁇ , which divides the period for which the advertisement is to be placed in units of 1 second, as the estimation target point 1000, and sets the estimation device 10 in place.
  • an estimated value 1500 that is, ⁇ (s 1 ), ⁇ (s 2 ), ..., ⁇ (s P )) of the purchase occurrence probability of the user at each time point in the period can be obtained.
  • the person in charge of advertising strategy can select users with a high probability of purchase occurrence on average under the hypothesis that users with a high probability of purchase occurrence have a high advertising effect, and can set up a schedule for placing advertisements with priority.
  • the parameter values may be determined by the parameter determination unit 103 described in the second embodiment.
  • the error calculation unit 104 described in the third embodiment may calculate the estimation error in addition to the purchase occurrence probability.
  • the rate function estimator 1400 It is output.
  • the person in charge sets the spatiotemporal point sequence ⁇ s 1 , s 2 , ..., S P ⁇ that divides the year into units of 1 second and the target area in units of 0.000001 degrees as the estimation target point 1000, and sets the estimation device.
  • an estimated value of 1500 crime rates at each spatiotemporal point that is, ⁇ (s 1 ), ⁇ (s 2 ), ..., ⁇ (s P )
  • the person in charge can formulate an effective annual allocation plan for security guards by identifying the time zone and area where the crime rate is high.
  • the infectious disease incidence rate it is possible to estimate the infectious disease incidence rate and utilize it by the same method as in the case of estimating the crime incidence rate for formulating an annual placement plan for security guards.
  • the past infectious disease outbreak history data of the target area may be used.
  • the infectious disease occurrence history data includes, for example, location information (latitude / longitude / altitude) of the target area, time information (time zone, year, etc.), weather information (sunny weather, rainy weather, etc.), and the like. (That is, the infectious disease occurrence history data may be, for example, a set of multidimensional vectors composed of information such as latitude, longitude, altitude, time zone, and weather). This makes it possible to estimate the rate function for the outbreak of infectious diseases under various conditions, and to formulate effective health and hygiene measures.
  • Non-Patent Document 1 As the conventional method to be compared with the proposed method, the method described in Non-Patent Document 1 above was adopted. Moreover, as an evaluation index, the average of the L2 norm errors (Average L2 norm errors) between the estimated value of the rate function and the true value was adopted. However, the average was the average of 10 trials. The evaluation results are shown in Table 1 below.
  • PP Gauss , PP Matern, and PP Wiener use the function h shown in the above equation (2) as the GCP kernel, and Gaussian kernel k Gauss (t, t') and Matern kernel k Matern (t, t'), respectively.
  • Wiener kernel k Wiener (t, t') is used as the semi-definite kernel k (t, t') to represent the proposed method.
  • VBPP represents the method described in Non-Patent Document 1 described above.
  • ⁇ 1 (t), ⁇ 2 (t) and ⁇ 3 (t) in Table 1 above are similar to those described in “4 Experiments” in the above reference and are included in the artificial dataset. This is the intensity function used to generate the point sequence.
  • the numbers in parentheses in Table 1 above are the standard errors of the L2 norm error between the estimated value and the true value of the rate function in 10 trials.
  • FIG. 7 is a diagram showing an example of the hardware configuration of the estimation device 10 according to the present embodiment.
  • the estimation device 10 is realized by a general computer or computer system, and includes an input device 201, a display device 202, an external I / F 203, a communication I / F 204, and a processor. It has 205 and a memory device 206. Each of these hardware is communicably connected via bus 207.
  • the input device 201 is, for example, a keyboard, a mouse, a touch panel, or the like.
  • the display device 202 is, for example, a display or the like.
  • the estimation device 10 does not have to have at least one of the input device 201 and the display device 202.
  • the external I / F 203 is an interface with an external device such as a recording medium 203a.
  • the estimation device 10 can read or write the recording medium 203a via the external I / F 203. Even if the recording medium 203a stores one or more programs that realize each functional unit (estimator generation unit 101, estimation unit 102, parameter determination unit 103, and error calculation unit 104) of the estimation device 10, for example. good.
  • the recording medium 203a includes, for example, a CD (Compact Disc), a DVD (Digital Versatile Disk), an SD memory card (Secure Digital memory card), a USB (Universal Serial Bus) memory card, and the like.
  • the communication I / F 204 is an interface for connecting the estimation device 10 to the communication network.
  • One or more programs that realize each functional unit of the estimation device 10 may be acquired (downloaded) from a predetermined server device or the like via the communication I / F 204.
  • the processor 205 is, for example, various arithmetic units such as a CPU (Central Processing Unit). Each functional unit included in the estimation device 10 is realized, for example, by a process in which one or more programs stored in the memory device 206 are executed by the processor 205.
  • CPU Central Processing Unit
  • the memory device 206 is, for example, various storage devices such as HDD, SSD, RAM (RandomAccessMemory), ROM (ReadOnlyMemory), and flash memory.
  • the estimation device 10 can realize the above-mentioned estimation process by having the hardware configuration shown in FIG. 7.
  • the hardware configuration shown in FIG. 7 is an example, and the estimation device 10 may have another hardware configuration.
  • the estimation device 10 may have a plurality of processors 205 or a plurality of memory devices 206.
  • the present invention is not limited to the above-described embodiment disclosed specifically, and various modifications and modifications, combinations with known techniques, and the like are possible without departing from the description of the claims. .. Moreover, the combination of each embodiment is also possible.
  • the estimation device 10 having the estimator generation unit 101, the estimation unit 102, the parameter determination unit 103, and the error calculation unit 104 may be configured by combining the second embodiment and the third embodiment.
  • Estimator 101 Estimator generator 102
  • Estimator 103 Parameter determination unit 104
  • Error calculation unit 201
  • Input device 202 Display device 203 External I / F 203a Recording medium 204 Communication I / F 205 Processor 206 Memory Device 207 Bus

Landscapes

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

Abstract

一実施形態に係る推定装置は、所定の次元の空間領域で観測されたイベントを表す観測点の集合と、前記観測が行われた回数を表す観測回数と、所定の条件を満たす第１の関数と、前記第１の関数のパラメータとを入力して、前記空間領域における前記イベントの発生率を得るためのレート関数を解析的に推定するレート関数推定部、を有することを特徴とする。

Description

推定装置、推定システム、推定方法及びプログラム
　本発明は、推定装置、推定システム、推定方法及びプログラムに関する。
　一次元又は多次元空間上で発生するイベントを表すデータを用いて、当該空間上の各点におけるイベントの発生率（以下、「レート関数」という。）を推定する技術として、Gaussian Cox Processが知られている。例えば、非特許文献１では、変分ベイズ法を用いた近似的なレート関数の推定方法が提案されている。
C. Lloyd, et al., "Variational Inference for Gaussian Process Modulated Poisson Processes", International Conference on Machine Learning, pp.1814-1822 (2015).
　しかしながら、変分ベイズ法を用いた近似は推定結果に誤差や偏りを生じさせる場合があり、真のレート関数と大きく異なる推定値が出力される可能性がある。また、その推定値の誤差の大きさを定量的に評価することが困難であり、誤った推定値を正しいものとして採用してしまう可能性もある。
　本発明の一実施形態は、上記の点に鑑みてなされたもので、解析的にレート関数を推定することを目的とする。
　上記目的を達成するため、一実施形態に係る推定装置は、所定の次元の空間領域で観測されたイベントを表す観測点の集合と、前記観測が行われた回数を表す観測回数と、所定の条件を満たす第１の関数と、前記第１の関数のパラメータとを入力して、前記空間領域における前記イベントの発生率を得るためのレート関数を解析的に推定するレート関数推定部、を有することを特徴とする。
　解析的にレート関数を推定することができる。
実施例１における推定装置の機能構成の一例を示す図である。 実施例１における推定処理の流れの一例を示すフローチャートである。 実施例２における推定装置の機能構成の一例を示す図である。 実施例２における推定処理の流れの一例を示すフローチャートである。 実施例３における推定装置の機能構成の一例を示す図である。 実施例３における推定処理の流れの一例を示すフローチャートである。 本実施形態に係る推定装置のハードウェア構成の一例を示す図である。
　以下、本発明の一実施形態について説明する。本実施形態では、Gaussian Cox Process（以下、「ＧＣＰ」と略記する。）に基づいて、解析的にレート関数を推定することができる推定装置１０について説明する。
　［実施例１］
　まず、本実施形態の実施例１について説明する。
　＜機能構成（実施例１）＞
　実施例１における推定装置１０の機能構成について、図１を参照しながら説明する。図１は、実施例１における推定装置１０の機能構成の一例を示す図である。
　図１に示すように、実施例１における推定装置１０は、推定器生成部１０１と、推定部１０２とを有する。また、実施例１における推定装置１０には、推定対象点１０００と、点事象データ１１００と、ＧＣＰカーネル１２００と、ＧＣＰカーネルパラメータ１３００とが与えられるものとする。なお、推定対象点１０００、点事象データ１１００、ＧＣＰカーネル１２００及びＧＣＰカーネルパラメータ１３００は、例えば、ＨＤＤ（Hard Disk Drive）やＳＳＤ（Solid State Drive）等の記憶装置に格納されていてもよいし、ユーザの操作等によって入力されてもよいし、通信ネットワークを介して接続される他の装置から受信してもよい。
　推定対象点１０００は、レート関数の値を推定する対象のＰ個の点（座標点）の集合｛ｓ１，ｓ２，・・・，ｓＰ｝であるものとする。これらの各ｓｐ（ｐ＝１，２，・・・，Ｐ）は、後述する空間領域Ｔ内で、レート関数の値の推定対象となる点のことである。
　点事象データ１１００は、観測を実施した空間領域Ｔと、その観測回数Ｒと、当該空間領域Ｔ内で観測されたＮ個のイベントの点（座標点）の集合｛ｔ１，ｔ２，・・・，ｔＮ｝である。なお、空間領域Ｔは予め決められた任意の次元（１次元又は多次元）の空間領域である。
　ＧＣＰカーネル１２００は、以下の条件１及び条件２を満たす関数ｈ（ｔ，ｔ'）である。ただし、ｔ，ｔ'はそれぞれ空間領域Ｔ内の任意の点である。
　条件１：ｈは対称な関数である。すなわち、ｈ（ｔ，ｔ'）＝ｈ（ｔ'，ｔ）である。
　条件２：関数ｈ（ｔ，ｔ'）を積分演算子のカーネルとして利用するとき、当該積分演算子の固有値が０以上１／Ｒ以下となる。ただし、関数ｈ（ｔ，ｔ'）の固有値νは以下の式（１）により定義される。
Figure JPOXMLDOC01-appb-M000001
 　ここで、φ（ｔ）は固有関数である。
　以降では、上記の条件１及び条件２を満たす関数を「ＧＣＰカーネル」又は「ＧＣＰカーネル関数」と呼ぶ。ＧＣＰカーネルの例を以下に示す。
　ＧＣＰカーネルの例１：半正定値カーネルｋ（ｔ，ｔ'）とその固有値の上限値ｋｍａｘで定義される以下の式（２）に示す関数ｈはＧＣＰカーネルとなる。
Figure JPOXMLDOC01-appb-M000002
 　ただし、ｋｍａｘは、以下の式（３）に示すように、半正定値カーネルｋ（ｔ，ｔ'）を一方の変数（式（３）ではｔ'）に関して積分した関数の最大値として評価される。
Figure JPOXMLDOC01-appb-M000003
 　半正定値カーネルｋ（ｔ，ｔ'）としては、強度αと共分散行列Σをパラメータに持つGaussianカーネルｋGauss（ｔ，ｔ'）や、強度αとスケースβをパラメータに持つMaternカーネルｋMatern（ｔ，ｔ'）等がある。GaussianカーネルｋGauss（ｔ，ｔ'）を以下の式（４）に、MaternカーネルｋMatern（ｔ，ｔ'）を以下の式（５）にそれぞれ示す。
Figure JPOXMLDOC01-appb-M000004
Figure JPOXMLDOC01-appb-M000005
 　ただし、τはベクトル及び行列に対する転置操作を表す。なお、GaussianカーネルｋGauss（ｔ，ｔ'）やMaternカーネルｋMatern（ｔ，ｔ'）以外にも、例えば、WienerカーネルｋWiener（ｔ，ｔ'）を用いることも可能である。
　ＧＣＰカーネルの例２：半正定値カーネルｋ（ｔ，ｔ'）を用いて定義される第二種フレドホルム型積分方程式の解はＧＣＰカーネルとなる。すなわち、以下の式（６）に示す積分方程式の解ｈ（ｔ，ｔ'）はＧＣＰカーネルとなる。
Figure JPOXMLDOC01-appb-M000006
 　ＧＣＰカーネルパラメータ１３００は、ＧＣＰカーネルｈ（ｔ，ｔ'）のパラメータ値である。ＧＣＰカーネルｈ（ｔ，ｔ'）のパラメータの個数及びパラメータが取り得る値の範囲は、ＧＣＰカーネルｈ（ｔ，ｔ'）の定義によって異なる。例えば、半正定値カーネルｋ（ｔ，ｔ'）にGaussianカーネルｋGauss（ｔ，ｔ'）を用いて、上記の式（２）に示す関数ｈをＧＣＰカーネルとした場合、そのパラメータは強度α及び共分散行列Σである。一方で、例えば、半正定値カーネルｋ（ｔ，ｔ'）にMaternカーネルｋMatern（ｔ，ｔ'）を用いて、上記の式（２）に示す関数ｈをＧＣＰカーネルとした場合、そのパラメータは強度α及びスケースβである。
　また、ＧＣＰカーネルパラメータ１３００には、レート関数の平均値の平方根を表すパラメータμが含まれる。パラメータμは非負の値を取るものとする。
　ＧＣＰカーネルパラメータ１３００に含まれる各パラメータの値は、通常、本実施形態に係る推定装置１０を適用する事例毎の様々な先行知識や制約等に基づいて決定される。例えば、ＥＣサイトのユーザが或る季節性の強い商材を購入するというイベントを対象として、MaternカーネルｋMatern（ｔ，ｔ'）を用いて、上記の式（２）に示す関数ｈをＧＣＰカーネルとした場合、パラメータβは数週程度を表す値とすることが考えられる。また、例えば、いくつかの候補となるパラメータ値を用いてレート関数をそれぞれ推定した上で、それらの推定結果を比較し、主観的に妥当と考えられるパラメータ値を採用することも考えられる。
　推定器生成部１０１は、推定装置１０に与えられた点事象データ１１００とＧＣＰカーネル１２００とＧＣＰカーネルパラメータ１３００とを入力して、レート関数推定器１４００を生成及び出力する。
　レート関数推定器１４００は、空間領域Ｔ内の任意の点ｔにおけるレート関数の推定値を得るための関数λ（ｔ）である。この関数λ（ｔ）は以下の式（７）で表される。
Figure JPOXMLDOC01-appb-M000007
 　ただし、
Figure JPOXMLDOC01-appb-M000008
 はＧＣＰカーネルｈ（ｔ，ｔ'）を空間領域Ｔにおいて積分した関数、ｈ（ｔ，：）はＮ個のイベントの座標点ｔ１，ｔ２，・・・，ｔＮ上におけるＧＣＰカーネルの値で構成されるベクトル、ｚはＮ個のイベントの座標点ｔ１，ｔ２，・・・，ｔＮ上におけるレート関数の推定値の平方根の逆数で構成されるベクトルとする。すなわち、
Figure JPOXMLDOC01-appb-M000009
 とする。
　上記の式（７）に示されるように、レート関数推定器１４００（つまり、λ（ｔ））は
Figure JPOXMLDOC01-appb-M000010
 の５つの要素で構成される。したがって、推定器生成部１０１は、これらの５つの要素を出力することで、レート関数推定器１４００を出力する。このとき、推定器生成部１０１は、これら５つの要素の中で未知数であるｚを計算により生成（算出）した上で、レート関数推定器１４００を出力する。ここで、Ｎ次元ベクトルｚの各要素ｚｎ＝λ－１／２（ｔｎ）はレート関数の推定値の関数であり、かつ、レート関数の推定値は上記の式（７）によりｚの関数であることから、この両者の関係が無矛盾に成り立つための条件式が存在し、その条件式を解くことでｚが算出される。この条件式は、以下の式（８）に示すＮ元連立二次方程式となる。
Figure JPOXMLDOC01-appb-M000011
 　なお、ｎ'も１≦ｎ'≦Ｎである。
　推定部１０２は、推定装置１０に与えられた推定対象点１０００と、推定器生成部１０１によって出力されたレート関数推定器１４００とを入力して、推定値１５００を算出及び出力する。すなわち、推定部１０２は、Ｐ個の点ｓ１，ｓ２，・・・，ｓＰのそれぞれを上記の式（７）に示すλ（ｔ）に代入することで、レート関数の推定値の集合｛λ（ｓ１），λ（ｓ２），・・・，λ（ｓＰ）｝を算出する。この｛λ（ｓ１），λ（ｓ２），・・・，λ（ｓＰ）｝が推定値１５００である。
　＜推定処理（実施例１）＞
　次に、実施例１における推定装置１０によって推定値１５００を得るための推定処理の流れについて、図２を参照しながら説明する。図２は、実施例１における推定処理の流れの一例を示すフローチャートである。
　ステップＳ１０１：まず、推定器生成部１０１は、与えられた点事象データ１１００とＧＣＰカーネル１２００とＧＣＰカーネルパラメータ１３００とを入力して、レート関数推定器１４００を生成及び出力する。
　ステップＳ１０２：そして、推定部１０２は、与えられた推定対象点１０００と、上記のステップＳ１０１で出力されたレート関数推定器１４００とを入力して、推定値１５００を算出及び出力する。これにより、レート関数の推定値の集合｛λ（ｓ１），λ（ｓ２），・・・，λ（ｓＰ）｝が得られる。このとき、実施例１における推定装置１０では、近似を用いずに、ＧＣＰに基づいて解析的にレート関数を推定するため、推定結果に意図しない誤差や偏りが含まれてしまう事態を防止することができる。
　なお、上記のステップＳ１０１の処理は、ステップＳ１０２の処理を実行する前に予め実行されていてもよい（この場合、レート関数推定器１４００を記憶装置等に保存しておき、上記のステップＳ１０１の処理では記憶装置等からレート関数推定器１４００を読み出せばよい。）。このとき、同一のレート関数推定器１４００を使用して上記のステップＳ１０２の処理を複数回実行する場合には、上記のステップＳ１０１の処理は１度実行されればよい。また、例えば、上記のステップＳ１０１の処理が他の装置で実行され、この他の装置からレート関数推定器１４００が与えられる場合には、推定装置１０は、上記のステップＳ１０２の処理のみを実行すればよい（この場合、当該他の装置と推定装置１０とでシステムが構成され、当該他の装置が推定器生成部１０１を有する一方で、推定装置１０は推定器生成部１０１を有していなくてもよい。）。
　［実施例２］
　次に、本実施形態の実施例２について説明する。実施例２では、ＧＣＰカーネルパラメータ１３００に含まれる各パラメータの値を主観的に決定するのではなく、主観を廃して決定する場合について説明する。
　＜機能構成（実施例２）＞
　実施例２における推定装置１０の機能構成について、図３を参照しながら説明する。図３は、実施例２における推定装置１０の機能構成の一例を示す図である。
　図３に示すように、実施例２における推定装置１０は、推定器生成部１０１と、推定部１０２と、パラメータ決定部１０３とを有する。また、実施例２における推定装置１０には、推定対象点１０００と、点事象データ１１００と、ＧＣＰカーネル１２００とが与えられるものとする。なお、推定対象点１０００、点事象データ１１００及びＧＣＰカーネル１２００は、例えば、ＨＤＤやＳＳＤ等の記憶装置に格納されていてもよいし、ユーザの操作等によって入力されてもよいし、通信ネットワークを介して接続される他の装置から受信してもよい。
　推定器生成部１０１及び推定部１０２は実施例１と同様であるため、その説明を省略する。ただし、推定器生成部１０１は、パラメータ決定部１０３によって出力されたＧＣＰカーネルパラメータ１３００を入力する。
　パラメータ決定部１０３は、推定装置１０に与えられた点事象データ１１００とＧＣＰカーネル１２００とを入力して、ＧＣＰカーネルパラメータ１３００を出力する。すなわち、パラメータ決定部１０３は、空間領域Ｔと観測回数ＲとＮ個のイベントの点の集合｛ｔ１，ｔ２，・・・，ｔＮ｝とＧＣＰカーネルｈ（ｔ，ｔ'）とを用いて、交差検定によりＧＣＰカーネルパラメータ１３００に含まれる各パラメータの値を最適化な値に決定する。具体的には、パラメータ決定部１０３は、例えば、Ｎ個のイベントの点の集合｛ｔ１，ｔ２，・・・，ｔＮ｝をＫ個のバッチ（部分集合）に分割し、Ｋ－１個のバッチを訓練用、残り１個のバッチを検証用として、検証用のバッチを入れ替えながらパラメータ値の変更と検証とを繰り返すことで、各パラメータ値を最適化する。なお、Ｋは空間領域Ｔや観測回数Ｒに応じて決定すればよい。
　＜推定処理（実施例２）＞
　次に、実施例２における推定装置１０によって推定値１５００を得るための推定処理の流れについて、図４を参照しながら説明する。図４は、実施例２における推定処理の流れの一例を示すフローチャートである。
　ステップＳ２０１：まず、パラメータ決定部１０３は、与えられた点事象データ１１００とＧＣＰカーネル１２００とを入力して、ＧＣＰカーネルパラメータ１３００を出力する。
　ステップＳ２０２：次に、推定器生成部１０１は、与えられた点事象データ１１００及びＧＣＰカーネル１２００と、上記のステップＳ２０１で出力されたＧＣＰカーネルパラメータ１３００とを入力して、レート関数推定器１４００を生成及び出力する。
　ステップＳ２０３：そして、推定部１０２は、与えられた推定対象点１０００と、上記のステップＳ２０２で出力されたレート関数推定器１４００とを入力して、推定値１５００を算出及び出力する。これにより、レート関数の推定値の集合｛λ（ｓ１），λ（ｓ２），・・・，λ（ｓＰ）｝が得られる。しかも、実施例２における推定装置１０では、ＧＣＰカーネルパラメータ１３００を最適化することができるため、例えば、本実施形態係る推定装置１０を適用する事例の先行知識や制約等に関する知識や経験が不足するユーザであっても、精度の良いレート関数の推定値を得ることが可能となる。
　なお、実施例１と同様に、上記のステップＳ２０１及びステップＳ２０２の処理が予め実行されていてもよいし、上記のステップＳ２０１及びステップＳ２０２の処理が他の装置で実行されてもよい。
　［実施例３］
　次に、本実施形態の実施例３について説明する。実施例３では、レート関数の推定値の誤差も算出する場合について説明する。
　＜機能構成（実施例３）＞
　実施例３における推定装置１０の機能構成について、図５を参照しながら説明する。図５は、実施例３における推定装置１０の機能構成の一例を示す図である。
　図５に示すように、実施例３における推定装置１０は、推定器生成部１０１と、推定部１０２と、誤差算出部１０４とを有する。また、実施例３における推定装置１０には、推定対象点１０００と、点事象データ１１００と、ＧＣＰカーネル１２００と、ＧＣＰカーネルパラメータ１３００とが与えられるものとする。なお、推定対象点１０００、点事象データ１１００、ＧＣＰカーネル１２００及びＧＣＰカーネルパラメータ１３００は、例えば、ＨＤＤやＳＳＤ等の記憶装置に格納されていてもよいし、ユーザの操作等によって入力されてもよいし、通信ネットワークを介して接続される他の装置から受信してもよい。
　推定器生成部１０１及び推定部１０２は実施例１と同様であるため、その説明を省略する。
　誤差算出部１０４は、推定装置１０に与えられた推定対象点１０００と、推定器生成部１０１によって出力されたレート関数推定器１４００と、推定部１０２によって出力された推定値１５００とを入力して、誤差１６００を算出及び出力する。誤差１６００は、レート関数の平方根に対し、推定対象とするＰ個の点ｓ１，ｓ２，・・・，ｓＰにおける推定誤差（共分散行列）である。これらの推定誤差σ（ｔ）は、以下の式（９）により算出される。
Figure JPOXMLDOC01-appb-M000012
 　ただし、Λ及びＨはそれぞれＮ×Ｎの行列であり、Λの（ｎ，ｎ'）成分をΛｎｎ'＝λ（ｔｎ）δｎｎ'、Ｈの（ｎ，ｎ'）成分をＨｎｎ'＝ｈ（ｔｎ，ｔｎ'）である。なお、δｎｎ'はｎ＝ｎ'のときは１、そうでないときは０を取る関数である。
　＜推定処理（実施例３）＞
　次に、実施例３における推定装置１０によって推定値１５００及び誤差１６００を得るための推定処理の流れについて、図６を参照しながら説明する。図６は、実施例３における推定処理の流れの一例を示すフローチャートである。
　ステップＳ３０１：まず、推定器生成部１０１は、与えられた点事象データ１１００とＧＣＰカーネル１２００とＧＣＰカーネルパラメータ１３００とを入力して、レート関数推定器１４００を生成及び出力する。
　ステップＳ３０２：次に、推定部１０２は、与えられた推定対象点１０００と、上記のステップＳ３０１で出力されたレート関数推定器１４００とを入力して、推定値１５００を算出及び出力する。
　ステップＳ３０３：そして、誤差算出部１０４は、与えられた推定対象点１０００と、上記のステップＳ３０１で出力されたレート関数推定器１４００と、上記のステップＳ３０２で出力された推定値１５００とを入力して、誤差１６００を算出及び出力する。これにより、レート関数の推定値の集合｛λ（ｓ１），λ（ｓ２），・・・，λ（ｓＰ）｝だけなく、推定誤差の集合｛σ（ｓ１），σ（ｓ２），・・・，σ（ｓＰ）｝も得られる。このため、実施例３における推定装置１０では、レート関数の推定値の誤差の大きさを定量的に評価することが可能となる。
　なお、実施例１と同様に、上記のステップＳ３０１の処理が予め実行されていてもよいし、上記のステップＳ３０１の処理が他の装置で実行されてもよい。
　＜応用例＞
　次に、本実施形態に係る推定装置１０を具体的な事例に適用した応用例について説明する。
　・応用例１
　ＥＣサイトにおける広告戦略担当者が、或る商材に関する次年度の夏季セール広告を優先的に送るユーザを特定する状況を考える。まず、ＥＣサイトの当該商材に関する購買履歴データ（ユーザＩＤと購買時刻とのペアで表されるデータの集合）からユーザ毎の購買時刻列を抽出した購買時刻データを準備する。１ユーザの購買時刻データが、イベントの点の集合｛ｔ１，ｔ２，・・・，ｔＮ｝となる。ここで、時刻は１秒単位の１年周期の値（01-01 00:00:00～12-31 23:59:59）を取るものとする。このとき、例えば、過去３年分の購買履歴データを用いるならば、空間領域Ｔは１次元空間、観測回数Ｒは３回となる。
　１ユーザの購買時刻データ、空間領域Ｔ及び観測回数Ｒが含まれる点事象データ１１００と、広告戦略担当者が選択したＧＣＰカーネル１２００及びＧＣＰカーネルパラメータ１３００とを推定装置１０に与えることで、レート関数推定器１４００が出力される。
　次に、広告戦略担当者は、広告を出す予定の期間を１秒単位で区切った時刻点列｛ｓ１，ｓ２，・・・，ｓＰ｝を推定対象点１０００として、推定装置１０に与える。これにより、当該期間の各時刻点における当該ユーザの購買発生確率の推定値１５００（つまり、λ（ｓ１），λ（ｓ２），・・・，λ（ｓＰ））が得られる。
　広告戦略担当者は、例えば、購買発生確率が高いユーザは広告効果が高いという仮説の下、購買発生確率が平均的に高いユーザを選別し、優先的に広告を出すスケジュールを組むことができる。
　なお、ＧＣＰカーネルのパラメータ値に関して業務上の先行知識がない場合は、実施例２で説明したパラメータ決定部１０３により、パラメータ値を決定してもよい。また、実施例３で説明した誤差算出部１０４により、購買発生確率に加えて、その推定誤差を算出してもよい。購買発生確率と共にその推定誤差を算出することで、例えば、各ユーザに対して、購買発生確率の推定値ではなく、購買発生確率が或る一定値以上を取る確率を算出することが可能となる。この確率は推定値の誤差も考慮するため、購買発生確率の高いユーザを特定するより頑強な指標として利用である可能性がある。
　・応用例２
　国や地方自治体あるいは警備会社が、対象地域における警備員の年間配置計画を策定する状況を考える。担当者は、まず、過去の犯罪発生履歴データ（発生時刻並びに発生位置の緯度及び経度で表される３次元ベクトルデータの集合）を準備する。この犯罪発生履歴データが、イベントの点の集合｛ｔ１，ｔ２，・・・，ｔＮ｝となる。ここで、時刻は１秒単位の１年周期の値（01-01 00:00:00～12-31 23:59:59）を取り、緯度及び経度はそれぞれ0.000001度単位の１０進数形式の値を取るものとする。このとき、例えば、過去５年分の犯罪発生履歴データを用いるならば、空間領域Ｔは３次元空間、観測回数Ｒは５回となる。
　犯罪発生履歴データ、空間領域Ｔ及び観測回数Ｒが含まれる点事象データ１１００と、担当者が選択したＧＣＰカーネル１２００及びＧＣＰカーネルパラメータ１３００とを推定装置１０に与えることで、レート関数推定器１４００が出力される。
　次に、担当者は、１年を１秒単位かつ対象地域を0.000001度単位で区切った時空間点列｛ｓ１，ｓ２，・・・，ｓＰ｝を推定対象点１０００として、推定装置１０に与える。これにより、各時空間点における犯罪発生率の推定値１５００（つまり、λ（ｓ１），λ（ｓ２），・・・，λ（ｓＰ））が得られる。したがって、担当者は、犯罪発生率が高い時間帯とエリアを特定することで、警備員の効果的な年間配置計画を策定することが可能となる。
　なお、警備員の年間配置計画策定のための犯罪発生率の推定事例と同様の方法により、例えば、感染症発生率の推定とその活用も可能である。この場合、対象地域の過去の感染症発生履歴データを用いればよい。このとき、感染症発生履歴データには、例えば、対象地域の位置情報（緯度・経度・高度）、時間情報（時間帯や年間等）、天候情報（晴天や雨天等）等が含まれていてもよい（つまり、感染症発生履歴データは、例えば、緯度、経度、高度、時間帯、天候等の情報で構成される多次元ベクトルの集合であってもよい。）。これにより、様々な条件下での感染症発生に関するレート関数を推定することができ、効果的な保健衛生に関する施策の策定が可能となる。
　＜評価＞
　次に、本実施形態に係る推定装置１０が実行する推定処理の手法（提案手法）の評価について説明する。本評価では、参考文献「Aglietti, V., Bonilla, E. V., Damoulas, T., and Cripps, S. Structured variational inference in continuous cox process models. In Advances in Neural Information Processing Systems, pp. 12437-12447, 2019.」の「4 Experiments」に記載されている方法と同様に、インターネット上の「ＵＲＬ：https://github.com/VirgiAgl/STVB」で提供されている人工データセットを用いて評価を行った。
　提案手法の比較対象とする従来手法としては、上記の非特許文献１に記載されている手法を採用した。また、評価指標としては、レート関数の推定値と真値とのＬ２ノルム誤差の平均（Average L2 norm errors）を採用した。ただし、平均は、１０回の試行の平均とした。以下の表１に評価結果を示す。
Figure JPOXMLDOC01-appb-T000013
 　ここで、ＰＰGauss、ＰＰMatern及びＰＰWienerは、上記の式（２）に示す関数ｈをＧＣＰカーネルとして、それぞれGaussianカーネルｋGauss（ｔ，ｔ'）、MaternカーネルｋMatern（ｔ，ｔ'）及びWienerカーネルｋWiener（ｔ，ｔ'）を半正定値カーネルｋ（ｔ，ｔ'）とした提案手法を表す。また、ＶＢＰＰは、上記の非特許文献１に記載されている手法を表す。
　上記の表１中のλ１（ｔ）、λ２（ｔ）及びλ３（ｔ）は上記の参考文献の「4 Experiments」に記載されているものと同様であり、人工データセットに含まれる点列の生成に用いられた強度関数である。なお、上記の表１中の括弧内は、１０回の試行における、レート関数の推定値と真値とのＬ２ノルム誤差の標準誤差（standard errors）である。
　上記の表１に示されるように、強度関数λ２（ｔ）で生成された点列を用いた場合のＰＰWienerはＶＢＰＰよりもやや精度が低いものの、他の場合では、提案手法は従来手法と比較して精度良くレート関数が推定できていることがわかる。
　＜ハードウェア構成＞
　最後に、本実施形態に係る推定装置１０のハードウェア構成について、図７を参照しながら説明する。図７は、本実施形態に係る推定装置１０のハードウェア構成の一例を示す図である。
　図７に示すように、本実施形態に係る推定装置１０は一般的なコンピュータ又はコンピュータシステムで実現され、入力装置２０１と、表示装置２０２と、外部Ｉ／Ｆ２０３と、通信Ｉ／Ｆ２０４と、プロセッサ２０５と、メモリ装置２０６とを有する。これら各ハードウェアは、それぞれがバス２０７を介して通信可能に接続されている。
　入力装置２０１は、例えば、キーボードやマウス、タッチパネル等である。表示装置２０２は、例えば、ディスプレイ等である。なお、推定装置１０は、入力装置２０１及び表示装置２０２のうちの少なくとも一方を有していなくてもよい。
　外部Ｉ／Ｆ２０３は、記録媒体２０３ａ等の外部装置とのインタフェースである。推定装置１０は、外部Ｉ／Ｆ２０３を介して、記録媒体２０３ａの読み取りや書き込み等を行うことができる。記録媒体２０３ａには、例えば、推定装置１０が有する各機能部（推定器生成部１０１、推定部１０２、パラメータ決定部１０３及び誤差算出部１０４）を実現する１以上のプログラムが格納されていてもよい。なお、記録媒体２０３ａには、例えば、ＣＤ（Compact Disc）、ＤＶＤ（Digital Versatile Disk）、ＳＤメモリカード（Secure Digital memory card）、ＵＳＢ（Universal Serial Bus）メモリカード等がある。
　通信Ｉ／Ｆ２０４は、推定装置１０を通信ネットワークに接続するためのインタフェースである。なお、推定装置１０が有する各機能部を実現する１以上のプログラムは、通信Ｉ／Ｆ２０４を介して、所定のサーバ装置等から取得（ダウンロード）されてもよい。
　プロセッサ２０５は、例えば、ＣＰＵ（Central Processing Unit）等の各種演算装置である。推定装置１０が有する各機能部は、例えば、メモリ装置２０６に格納されている１以上のプログラムがプロセッサ２０５に実行させる処理により実現される。
　メモリ装置２０６は、例えば、ＨＤＤやＳＳＤ、ＲＡＭ（Random Access Memory）、ＲＯＭ（Read Only Memory）、フラッシュメモリ等の各種記憶装置である。
　本実施形態に係る推定装置１０は、図７に示すハードウェア構成を有することにより、上述した推定処理を実現することができる。なお、図７に示すハードウェア構成は一例であって、推定装置１０は、他のハードウェア構成を有していてもよい。例えば、推定装置１０は、複数のプロセッサ２０５を有していてもよいし、複数のメモリ装置２０６を有していてもよい。
　本発明は、具体的に開示された上記の実施形態に限定されるものではなく、請求の範囲の記載から逸脱することなく、種々の変形や変更、既知の技術との組み合わせ等が可能である。また、各実施例の組み合わせも可能である。例えば、実施例２と実施例３とを組み合わせて、推定器生成部１０１と推定部１０２とパラメータ決定部１０３と誤差算出部１０４とを有する推定装置１０が構成されてもよい。
　１０　　　　推定装置
　１０１　　　推定器生成部
　１０２　　　推定部
　１０３　　　パラメータ決定部
　１０４　　　誤差算出部
　２０１　　　入力装置
　２０２　　　表示装置
　２０３　　　外部Ｉ／Ｆ
　２０３ａ　　記録媒体
　２０４　　　通信Ｉ／Ｆ
　２０５　　　プロセッサ
　２０６　　　メモリ装置
　２０７　　　バス

Claims (7)

  1. 　所定の次元の空間領域で観測されたイベントを表す観測点の集合と、前記観測が行われた回数を表す観測回数と、所定の条件を満たす第１の関数と、前記第１の関数のパラメータとを入力して、前記空間領域における前記イベントの発生率を得るためのレート関数を解析的に推定するレート関数推定部、
    　を有することを特徴とする推定装置。
  2. 　前記所定の条件は、
    　前記第１の関数が対称関数であること、かつ、前記観測回数をＲとして、前記第１の関数を積分演算子のカーネル関数として利用した場合における固有値が０以上１／Ｒ以下であること、を特徴とする請求項１に記載の推定装置。
  3. 　前記空間領域内の点であり、かつ、前記イベントの発生率の推定対象となる点を表す推定対象点の集合と、前記レート関数とを入力して、前記推定対象点におけるイベントの発生率を推定する発生率推定部、を更に有することを特徴とする請求項１又は２に記載の推定装置。
  4. 　前記レート関数推定部は、
    　前記観測回数と、前記レート関数の平均値の平方根を表すパラメータと、前記第１の関数を前記空間領域で積分した第２の関数と、前記観測点における前記第１の関数の値を表すベクトルと、前記観測点における前記レート関数の推定値の平方根の逆数を表すベクトルとが含まれる前記レート関数を解析的に推定する、ことを特徴とする請求項１乃至３の何れか一項に記載の推定装置。
  5. 　第１の装置と、前記第１の装置とは異なる第２の装置とが含まれる推定システムであって、
    　前記第１の装置は、
    　所定の次元の空間領域で観測されたイベントを表す観測点の集合と、前記観測が行われた回数を表す観測回数と、所定の条件を満たす第１の関数と、前記第１の関数のパラメータとを入力して、前記空間領域における前記イベントの発生率を得るためのレート関数を解析的に推定するレート関数推定部、を有し、
    　前記第２の装置は、
    　前記空間領域内の点であり、かつ、前記イベントの発生率の推定対象となる点を表す推定対象点の集合と、前記レート関数とを入力して、前記推定対象点におけるイベントの発生率を推定する発生率推定部、を有することを特徴とする推定システム。
  6. 　所定の次元の空間領域で観測されたイベントを表す観測点の集合と、前記観測が行われた回数を表す観測回数と、所定の条件を満たす第１の関数と、前記第１の関数のパラメータとを入力して、前記空間領域における前記イベントの発生率を得るためのレート関数を解析的に推定するレート関数推定手順、
    　をコンピュータが実行することを特徴とする推定方法。
  7. 　コンピュータを、請求項１乃至４の何れか一項に記載の推定装置、又は、請求項５に記載の推定システムに含まれる第１の装置若しくは第２の装置、として機能させるプログラム。
PCT/JP2020/018115 2020-04-28 2020-04-28 推定装置、推定システム、推定方法及びプログラム WO2021220399A1 (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
US17/919,910 US20230205841A1 (en) 2020-04-28 2020-04-28 Estimation apparatus, estimation system, estimation method and program
JP2022518485A JP7392838B2 (ja) 2020-04-28 2020-04-28 推定装置、推定システム、推定方法及びプログラム
PCT/JP2020/018115 WO2021220399A1 (ja) 2020-04-28 2020-04-28 推定装置、推定システム、推定方法及びプログラム

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
PCT/JP2020/018115 WO2021220399A1 (ja) 2020-04-28 2020-04-28 推定装置、推定システム、推定方法及びプログラム

Publications (1)

Publication Number Publication Date
WO2021220399A1 true WO2021220399A1 (ja) 2021-11-04

Family

ID=78373459

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2020/018115 WO2021220399A1 (ja) 2020-04-28 2020-04-28 推定装置、推定システム、推定方法及びプログラム

Country Status (3)

Country Link
US (1) US20230205841A1 (ja)
JP (1) JP7392838B2 (ja)
WO (1) WO2021220399A1 (ja)

Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2019008637A (ja) * 2017-06-27 2019-01-17 日本電信電話株式会社 パラメータ推定装置、機器故障予測装置、最適施策探索装置、方法、及びプログラム

Patent Citations (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2019008637A (ja) * 2017-06-27 2019-01-17 日本電信電話株式会社 パラメータ推定装置、機器故障予測装置、最適施策探索装置、方法、及びプログラム

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
SAWAMURA, YASUMASA; AKIRA, HAYASHI; NOBUO, SUEMATSU: "Fast similarity search of time series data using the nystrom method", 5TH FORUM ON INFORMATION TECHNOLOGY (FIT2006) GENERAL LECTURE PROCEEDINGS; FUKUOKA; SEPTEMBER 5-7, 2006, vol. 5, no. 3, 21 August 2006 (2006-08-21), JP, pages 51 - 54, XP009532011 *

Also Published As

Publication number Publication date
JP7392838B2 (ja) 2023-12-06
JPWO2021220399A1 (ja) 2021-11-04
US20230205841A1 (en) 2023-06-29

Similar Documents

Publication Publication Date Title
CN110427560B (zh) 一种应用于推荐系统的模型训练方法以及相关装置
JP6414363B2 (ja) 予測システム、方法およびプログラム
US10515378B2 (en) Extracting relevant features from electronic marketing data for training analytical models
Wang et al. Cost-effective quality assurance in crowd labeling
Keane et al. A computationally practical simulation estimation algorithm for dynamic panel data models with unobserved endogenous state variables
TW201734893A (zh) 信用分的獲取、特徵向量值的輸出方法及其裝置
US10984343B2 (en) Training and estimation of selection behavior of target
JP6311851B2 (ja) 共クラスタリングシステム、方法およびプログラム
CN113407987B (zh) 保护隐私的确定业务数据特征有效值的方法及装置
KR102192949B1 (ko) 인공신경망을 이용한 창업 기업의 분석 방법 및 장치
WO2014103560A1 (ja) 分析装置、分析プログラム、分析方法、推定装置、推定プログラム、及び、推定方法。
US20210151128A1 (en) Learning Method, Mixing Ratio Prediction Method, and Prediction Device
US20230401426A1 (en) Prediction method, prediction apparatus and program
Zhou et al. Multiple imputation in two-stage cluster samples using the weighted finite population Bayesian bootstrap
Zhu et al. Weak proxies are sufficient and preferable for fairness with missing sensitive attributes
Koduvely Learning Bayesian Models with R
Tong et al. Learning fractional white noises in neural stochastic differential equations
WO2021220399A1 (ja) 推定装置、推定システム、推定方法及びプログラム
Sharipov et al. Bootstrap for the sample mean and for U-statistics of mixing and near-epoch dependent processes
De Neve et al. Goodness-of-fit methods for probabilistic index models
Vinzamuri et al. Pre-processing censored survival data using inverse covariance matrix based calibration
Kravdal The problematic estimation of “imitation effects” in multilevel models
Mikhail et al. Estimating persistence in Canadian unemployment: evidence from a Bayesian ARFIMA
Mbaluka et al. Application of principal component analysis and hierarchical regression model on Kenya macroeconomic indicators
Kürüm et al. Multilevel joint modeling of hospitalization and survival in patients on dialysis

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

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2022518485

Country of ref document: JP

Kind code of ref document: A

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 20933907

Country of ref document: EP

Kind code of ref document: A1