WO2020026327A1 - 情報処理装置、制御方法、及びプログラム - Google Patents

情報処理装置、制御方法、及びプログラム Download PDF

Info

Publication number
WO2020026327A1
WO2020026327A1 PCT/JP2018/028565 JP2018028565W WO2020026327A1 WO 2020026327 A1 WO2020026327 A1 WO 2020026327A1 JP 2018028565 W JP2018028565 W JP 2018028565W WO 2020026327 A1 WO2020026327 A1 WO 2020026327A1
Authority
WO
WIPO (PCT)
Prior art keywords
time
series data
feature
contribution
sensor
Prior art date
Application number
PCT/JP2018/028565
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 JP2020533925A priority Critical patent/JP7140191B2/ja
Priority to US17/262,955 priority patent/US20210311009A1/en
Priority to PCT/JP2018/028565 priority patent/WO2020026327A1/ja
Publication of WO2020026327A1 publication Critical patent/WO2020026327A1/ja

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N33/00Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
    • G01N33/0004Gaseous mixtures, e.g. polluted air
    • G01N33/0009General constructional details of gas analysers, e.g. portable test equipment
    • G01N33/0062General constructional details of gas analysers, e.g. portable test equipment concerning the measuring method, e.g. intermittent, or the display, e.g. digital
    • G01N33/0067General constructional details of gas analysers, e.g. portable test equipment concerning the measuring method, e.g. intermittent, or the display, e.g. digital by measuring the rate of variation of the concentration
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N1/00Sampling; Preparing specimens for investigation
    • G01N1/02Devices for withdrawing samples
    • G01N1/22Devices for withdrawing samples in the gaseous state
    • G01N1/24Suction devices
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N19/00Investigating materials by mechanical methods
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N33/00Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
    • G01N33/0004Gaseous mixtures, e.g. polluted air
    • G01N33/0009General constructional details of gas analysers, e.g. portable test equipment
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N33/00Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
    • G01N33/0004Gaseous mixtures, e.g. polluted air
    • G01N33/0009General constructional details of gas analysers, e.g. portable test equipment
    • G01N33/0011Sample conditioning
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N33/00Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
    • G01N33/0004Gaseous mixtures, e.g. polluted air
    • G01N33/0009General constructional details of gas analysers, e.g. portable test equipment
    • G01N33/0011Sample conditioning
    • G01N33/0016Sample conditioning by regulating a physical variable, e.g. pressure, temperature
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N33/00Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
    • G01N33/0004Gaseous mixtures, e.g. polluted air
    • G01N33/0009General constructional details of gas analysers, e.g. portable test equipment
    • G01N33/0027General constructional details of gas analysers, e.g. portable test equipment concerning the detector
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N33/00Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
    • G01N33/0004Gaseous mixtures, e.g. polluted air
    • G01N33/0009General constructional details of gas analysers, e.g. portable test equipment
    • G01N33/0027General constructional details of gas analysers, e.g. portable test equipment concerning the detector
    • G01N33/0031General constructional details of gas analysers, e.g. portable test equipment concerning the detector comprising two or more sensors, e.g. a sensor array
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N5/00Analysing materials by weighing, e.g. weighing small particles separated from a gas or liquid
    • G01N5/02Analysing materials by weighing, e.g. weighing small particles separated from a gas or liquid by absorbing or adsorbing components of a material and determining change of weight of the adsorbent, e.g. determining moisture content

Definitions

  • the present invention relates to the analysis of gas characteristics.
  • Patent Literature 1 discloses a technique for determining the type of a sample gas using a signal (time-series data of a detection value) obtained by measuring the sample gas with a nanomechanical sensor. Specifically, the diffusion time constant of the sample gas with respect to the receptor of the sensor is determined by the combination of the type of the receptor and the type of the sample gas, and therefore, based on the diffusion time constant obtained from the signal and the type of the receptor. It is disclosed that the type of the sample gas can be determined.
  • Patent Document 1 assumes that the sample gas contains only one type of molecule, and does not assume that a sample gas in which a plurality of types of molecules are mixed is handled.
  • the present invention has been made in view of the above problems, and has as its object to provide a technique for extracting characteristics of a gas in which a plurality of types of molecules are mixed.
  • the information processing apparatus includes: 1) a time-series data acquisition unit configured to acquire time-series data of a detection value, which is output from a sensor whose detection value changes according to attachment and detachment of molecules contained in a target gas; 2) a calculation unit that calculates a plurality of feature constants that contribute to the time-series data and a contribution value indicating the magnitude of the contribution of each feature constant to the time-series data; 3) a plurality of feature constants and each of the feature constants And an output unit that outputs a combination of the contribution values calculated with respect to as a feature amount of the gas sensed by the sensor.
  • the characteristic constant is a time constant or a rate constant relating to the magnitude of a temporal change in the amount of molecules attached to the sensor.
  • the control method of the present invention is executed by a computer.
  • the control method includes: 1) a time-series data acquisition step of acquiring time-series data of a detection value output from a sensor whose detection value changes in accordance with the attachment and detachment of molecules contained in the target gas; A calculating step of calculating a plurality of feature constants contributing to the series data and a contribution value indicating the magnitude of the contribution of each feature constant to the time-series data; 3) calculating a plurality of feature constants and each feature constant An output step of outputting the combination of the contribution values obtained as the characteristic amount of the gas sensed by the sensor.
  • the characteristic constant is a time constant or a rate constant relating to the magnitude of a temporal change in the amount of molecules attached to the sensor.
  • the program of the present invention causes a computer to execute each step of the control method of the present invention.
  • FIG. 2 is a diagram illustrating an outline of the information processing apparatus according to the first embodiment.
  • FIG. 4 is a diagram illustrating a sensor for obtaining data acquired by the information processing device.
  • FIG. 2 is a diagram illustrating a functional configuration of the information processing apparatus according to the first embodiment.
  • FIG. 2 is a diagram illustrating a computer for realizing an information processing device.
  • 6 is a flowchart illustrating a flow of a process executed by the information processing apparatus according to the first embodiment.
  • FIG. 3 is a diagram illustrating a plurality of time-series data obtained from a sensor. It is a figure which illustrates the graph of g (t).
  • FIG. 7 is a diagram illustrating a case where a feature matrix is obtained from each of rising time series data and falling time series data. It is a figure which illustrates the case where a plurality of feature matrices are obtained by obtaining time-series data from each of a plurality of sensors.
  • each block diagram represents a configuration of a functional unit, not a configuration of a hardware unit.
  • FIG. 1 is a diagram illustrating an outline of an information processing apparatus 2000 according to the first embodiment.
  • FIG. 2 is a diagram illustrating a sensor 10 for obtaining data acquired by the information processing apparatus 2000.
  • the sensor 10 has a receptor to which a molecule is attached, and a detection value changes according to attachment and detachment of the molecule at the receptor.
  • the gas sensed by the sensor 10 is referred to as a target gas.
  • the time-series data of the detection values output from the sensor 10 is referred to as time-series data 14.
  • the time-series data 14 is also described as Y
  • the detected value at the time t is also described as y (t).
  • Y is a vector in which y (t) is enumerated.
  • the senor 10 is a membrane-type surface stress (MSS) sensor.
  • the MSS sensor has, as a receptor, a functional film to which a molecule adheres, and the stress generated in a support member of the functional film changes due to the attachment and detachment of the molecule to and from the functional film.
  • the MSS sensor outputs a detection value based on the change in the stress.
  • the sensor 10 is not limited to the MSS ⁇ sensor, and the physical quantity related to the viscoelasticity and dynamic characteristics (mass, moment of inertia, etc.) of the members of the sensor 10 that occur in response to the attachment and detachment of the molecule to and from the receptor. Any sensor that outputs a detection value based on the change may be used, and various types of sensors such as a cantilever type, a film type, an optical type, a piezo, and a vibration response can be employed.
  • the sensing by the sensor 10 is modeled as follows. (1) The sensor 10 is exposed to a target gas containing K kinds of molecules. (2) The concentration of each molecule k in the target gas is constant ⁇ k. (3) The sensor 10 can adsorb a total of N molecules. (4) The number of molecules k attached to the sensor 10 at time t is nk (t).
  • the time change of the number nk (t) of the molecules k attached to the sensor 10 can be formulated as follows.
  • the first and second terms on the right-hand side of the equation (1) are the increasing amount of the molecule ⁇ k ⁇ per unit time (the number of molecules ⁇ k ⁇ newly attached to the sensor 10) and the decreasing amount (the molecule ⁇ k ⁇ detached from the sensor 10). Number).
  • ⁇ k and ⁇ k are a rate constant representing the rate at which the molecule ⁇ k ⁇ adheres to the sensor 10 and a rate constant representing the rate at which the molecule ⁇ k ⁇ separates from the sensor 10, respectively.
  • the concentration ⁇ k is constant
  • the number nk (t) of the numerator k at time t can be formulated from the above equation (1) as follows.
  • nk (t) is expressed as follows.
  • the detection value of the sensor 10 is determined by the stress applied to the sensor 10 by molecules contained in the target gas. Then, it is considered that the stress acting on the sensor 10 by a plurality of molecules can be represented by a linear sum of the stress acting on each molecule. However, it is considered that the stress generated by the molecule differs depending on the type of the molecule. That is, it can be said that the contribution of the molecule to the detection value of the sensor 10 differs depending on the type of the molecule.
  • the detection value y (t) of the sensor 10 can be formulated as follows.
  • both ⁇ k and ⁇ k represent the contribution of the numerator k to the detection value of the sensor 10. The meaning of “rising” and “falling” will be described later.
  • the time-series data 14 obtained from the sensor 10 sensing the target gas can be decomposed as in the above equation (4), the types of molecules contained in the target gas and each type of molecule are contained in the target gas.
  • the information processing device 2000 acquires the time-series data 14 output by the sensor 10, and decomposes the time-series data 14 as shown in the following equation (5).
  • ⁇ i is a contribution value representing the contribution of the characteristic constant ⁇ i to the detection value of the sensor 10.
  • the contribution value ⁇ i ⁇ representing the magnitude of the contribution of the constant ⁇ i to the time-series data 14 is calculated.
  • the contribution value ⁇ i ⁇ is calculated after calculating the characteristic constant ⁇ i, and a case where the contribution value ⁇ i ⁇ is calculated together with the characteristic constant ⁇ i.
  • the information processing apparatus 2000 outputs information in which the set of feature constants and the set of contribution values are associated with each other as a feature amount representing the feature of the target gas.
  • the association between the set of feature constants ⁇ ⁇ and the set of contribution values ⁇ is represented, for example, by m feature matrix ⁇ F ⁇ having two rows and two columns (m is the number of each of the feature constant and the contribution value).
  • the characteristic amount of the target gas is represented by a characteristic matrix ⁇ F ⁇ unless otherwise specified.
  • the feature quantity of the target gas does not necessarily have to be represented as a vector.
  • equation (5) can be expressed as follows.
  • the velocity constant ⁇ is used as the characteristic constant. Therefore, the set ⁇ ⁇ of feature constants indicates ⁇ 1 to ⁇ m.
  • the above-described set of feature constants ⁇ and the set of corresponding contribution values ⁇ are included in the target gas. It is thought that they differ depending on the type of molecules and their mixing ratio. Therefore, the information in which the set of feature constants ⁇ and the set of contribution values ⁇ ⁇ are associated with each other can be used as information capable of distinguishing a gas in which a plurality of types of molecules are mixed, that is, a gas feature amount. it can.
  • the information processing apparatus 2000 of the present embodiment calculates a set ⁇ of characteristic constants and the contribution of each characteristic constant to the time-series data 14 based on the time-series data 14 obtained by sensing the target gas with the sensor 10.
  • a set ⁇ ⁇ of contribution values to be represented is calculated, and information in which the calculated sets ⁇ and ⁇ are associated with each other is output as a feature amount of the target gas. This makes it possible to automatically generate a characteristic amount capable of identifying a gas in which a plurality of types of molecules are mixed from the result of sensing the gas with the sensor 10.
  • FIG. 1 The above description with reference to FIG. 1 is an example for facilitating understanding of the information processing device 2000, and does not limit the functions of the information processing device 2000.
  • the information processing apparatus 2000 of the present embodiment will be described in more detail.
  • FIG. 3 is a diagram illustrating a functional configuration of the information processing apparatus 2000 according to the first embodiment.
  • the information processing device 2000 includes a time-series data acquisition unit 2020, a calculation unit 2040, and an output unit 2060.
  • the time-series data acquisition unit 2020 acquires the time-series data 14 from the sensor 10.
  • the feature constant generation unit 2030 uses the time series data 14 to calculate a plurality of feature constants ⁇ i and a contribution value ⁇ i representing the magnitude of the contribution of each feature constant ⁇ i to the time series data 14.
  • the output unit 2060 outputs information associating the set of feature constants ⁇ with the set of contribution values ⁇ ⁇ as the feature amount of the gas sensed by the sensor 10.
  • Each functional component of the information processing apparatus 2000 may be implemented by hardware (eg, a hard-wired electronic circuit or the like) that implements each functional component, or a combination of hardware and software (eg: Electronic circuit and a program for controlling the same).
  • hardware eg, a hard-wired electronic circuit or the like
  • software eg: Electronic circuit and a program for controlling the same.
  • FIG. 4 is a diagram illustrating a computer 1000 for realizing the information processing device 2000.
  • the computer 1000 is an arbitrary computer.
  • the computer 1000 is a stationary computer such as a personal computer (PC) or a server machine.
  • the computer 1000 is a portable computer such as a smartphone or a tablet terminal.
  • the computer 1000 may be a dedicated computer designed to realize the information processing device 2000, or may be a general-purpose computer.
  • the computer 1000 has a bus 1020, a processor 1040, a memory 1060, a storage device 1080, an input / output interface 1100, and a network interface 1120.
  • the bus 1020 is a data transmission path through which the processor 1040, the memory 1060, the storage device 1080, the input / output interface 1100, and the network interface 1120 mutually transmit and receive data.
  • a method for connecting the processors 1040 and the like to each other is not limited to a bus connection.
  • the processor 1040 is various processors such as a CPU (Central Processing Unit), a GPU (Graphics Processing Unit), and an FPGA (Field-Programmable Gate Array).
  • the memory 1060 is a main storage device realized using a RAM (Random Access Memory) or the like.
  • the storage device 1080 is an auxiliary storage device realized using a hard disk, an SSD (Solid State Drive), a memory card, or a ROM (Read Only Memory).
  • the input / output interface 1100 is an interface for connecting the computer 1000 and an input / output device.
  • an input device such as a keyboard and an output device such as a display device are connected to the input / output interface 1100.
  • the sensor 10 is connected to the input / output interface 1100.
  • the sensor 10 does not necessarily need to be directly connected to the computer 1000.
  • the sensor 10 may store the time-series data 14 in a storage device shared with the computer 1000.
  • the network interface 1120 is an interface for connecting the computer 1000 to a communication network.
  • the communication network is, for example, a LAN (Local Area Network) or a WAN (Wide Area Network).
  • the method by which the network interface 1120 connects to the communication network may be a wireless connection or a wired connection.
  • the storage device 1080 stores a program module that implements each functional component of the information processing apparatus 2000.
  • the processor 1040 realizes a function corresponding to each program module by reading out each of these program modules into the memory 1060 and executing them.
  • FIG. 5 is a flowchart illustrating a flow of a process executed by the information processing apparatus 2000 according to the first embodiment.
  • the time-series data acquisition unit 2020 acquires the time-series data 14 (S102).
  • the calculating unit 2040 calculates a plurality of characteristic constants ⁇ i and a contribution value ⁇ i corresponding to each characteristic constant ⁇ i using the time series data 14 (S104).
  • the output unit 2060 outputs information in which the set of feature constants and the set of contribution values are associated with each other as the feature amount of the target gas (S106).
  • the information processing apparatus 2000 executes the series of processes shown in FIG.
  • the information processing device 2000 receives an input operation for specifying the time-series data 14 and executes a series of processes on the specified time-series data 14.
  • the information processing apparatus 2000 stands by so as to be able to receive the time-series data 14, and in response to receiving the time-series data 14 (that is, executing S102), performs the processing after S104.
  • the time-series data acquisition unit 2020 acquires the time-series data 14 (S102).
  • the method by which the time-series data acquisition unit 2020 acquires the time-series data 14 is arbitrary.
  • the information processing apparatus 2000 acquires the time-series data 14 by accessing a storage device in which the time-series data 14 is stored.
  • the storage device in which the time series data 14 is stored may be provided inside the sensor 10 or may be provided outside the sensor 10.
  • the time-series data acquisition unit 2020 may obtain the time-series data 14 by sequentially receiving the detection values output from the sensor 10.
  • the time-series data 14 is time-series data in which the detection values output by the sensor 10 are arranged in ascending order of the time output from the sensor 10.
  • the time-series data 14 may be obtained by subjecting the time-series data of the detection value obtained from the sensor 10 to predetermined preprocessing.
  • the time-series data acquisition unit 2020 may perform pre-processing on the time-series data 14.
  • the pre-processing for example, filtering for removing a noise component from the time-series data can be employed.
  • the time-series data 14 is obtained by exposing the sensor 10 to the target gas.
  • the operation of exposing the sensor to the gas to be measured and the operation of removing the gas to be measured from the sensor are repeated to obtain a plurality of time-series data to be analyzed from the sensor. There is.
  • FIG. 6 is a diagram illustrating a plurality of time-series data obtained from sensors.
  • the rising time series data is represented by a solid line
  • the falling time series data is represented by a dotted line so that the rising time series data and the falling time series data can be easily distinguished.
  • the time series data 14-1 in the period P1 and the time series data 14-3 in the period P3 are obtained by an operation of exposing the sensor to the gas to be measured.
  • time series data obtained by exposing the sensor to the gas to be measured is also referred to as “rise” time series data.
  • “In the case of rising” in Expression (4) means “in the case where the time-series data 14 is time-series data of rising”. The same applies to the following equations.
  • the time series data 14-2 in the period P2 and the time series data 14-4 in the period P4 are obtained by an operation of removing the gas to be measured from the sensor.
  • the operation of removing the gas to be measured from the sensor is realized, for example, by exposing the sensor to a gas called a purge gas.
  • the time-series data obtained by the operation of removing the gas to be measured from the sensor is also referred to as “falling” time-series data.
  • “In the case of falling” in Equation (4) means “in the case where the time-series data 14 is falling time-series data”. The same applies to the following equations.
  • the time series data 14 obtained by the operation of exposing the sensor 10 to the target gas and the operation of removing the target gas from the sensor 10 are distinguished and treated as different time series data 14.
  • the time series data obtained in each of the four periods P1 to P4 is treated as different time series data 14. Therefore, when a series of time-series data is obtained by repeating the operation of exposing the sensor 10 to the target gas and the operation of removing the target gas from the sensor 10, the series of time-series data is converted into a plurality of time-series data 14. Need to be split.
  • a plurality of time series data 14 is obtained by manually dividing a series of time series data obtained from the sensor 10.
  • the information processing apparatus 2000 may acquire a series of time-series data, and divide the time-series data to obtain a plurality of time-series data 14.
  • ⁇ (c) Method of using metadata obtained from sensor
  • metadata other than the detected value is provided.
  • different pumps sample pump and purge pump
  • the recorded detection values include the operation sequence of the pump (information indicating which pump is used to detect the value, a flow rate measurement value used for feedback control of the flow rate, etc.) in time series. It is added as information. Therefore, for example, the information processing apparatus 2000 can divide the time-series data 14 by using the operation sequence of the pump obtained together with the time-series data 14.
  • the information processing apparatus 2000 tentatively divides the time-series data 14 into a plurality of sections by the method (c), specifies a time point at which the absolute value of the first-order differentiation becomes maximum in each section, and specifies The time series data 14 is divided at each time point.
  • the information processing apparatus 2000 uses only one of the time-series data 14 obtained by the operation of exposing the sensor 10 to the target gas and the time-series data 14 obtained by the operation of removing the target gas from the sensor 10. It may be configured as follows.
  • the calculation unit 2040 calculates a plurality of characteristic constants ⁇ i and a contribution value ⁇ i corresponding to the characteristic constant ⁇ i using the time-series data 14 (S104). This processing is equivalent to decomposing the time-series data 14 into the sum of ⁇ i * f ( ⁇ i) shown in Expression (5).
  • Various methods can be used as a method for calculating the feature constant and its contribution value from the time-series data 14 by the calculation unit 2040.
  • a method of calculating the feature constant will be described first.
  • the contribution value is calculated in the process of calculating the characteristic constant. This case will be described later.
  • Characteristic constant calculation method 1 For example, the calculation unit 2040 calculates a set of feature constants based on the slope of log [y '(t)] obtained by taking log as the differential y' (t) of the time-series data 14. log [y '(t)] is as follows.
  • the calculation unit 2040 calculates the gradient ⁇ g '(t) ⁇ of ⁇ g (t) ⁇ at each time ⁇ t ⁇ , and from the domain of g (t) (the measurement period of the time-series data 14), g ′ (t) is approximately A plurality of identical periods (hereinafter, partial periods) are extracted. Then, the calculating unit 2040 calculates ⁇ g ′ (t) ⁇ in each of the extracted partial periods as a speed constant ⁇ i ⁇ . For example, the calculating unit 2040 calculates a statistical value (such as an average value) of ⁇ g ′ (t) ⁇ in each partial period as a speed constant corresponding to the partial period.
  • the calculating unit 2040 performs a linear regression on ⁇ (t, ⁇ g (t)) ⁇ included in the partial period for each partial period, and calculates a slope of a regression line obtained for each partial period as a rate constant. Is also good. Note that a period in which g ′ (t) is substantially the same can be specified by, for example, clustering a set of g ′ (t) with the magnitude of the value.
  • a gradient different from any of the gradients before and after the gradient may appear in a portion where the gradient of g (t) is switched for a short period of time.
  • a slope that appears only for such a short period may not be included in the feature constant.
  • the minimum value of the length of the partial period is determined in advance. Then, the calculating unit 2040 divides the domain of ⁇ g (t) ⁇ into a plurality of periods based on the magnitude of the gradient ⁇ g ′ (t) ⁇ , and sets only the period having a length equal to or longer than the minimum value to the aforementioned partial period (Ie, ⁇ g ′ (t) ⁇ corresponding to a period having a length less than the minimum value is not included in the set of feature constants). By doing so, ⁇ g '(t) ⁇ that appears only for a short period can be excluded from the feature constant. Note that the minimum value of the length of the partial period can be expressed as the minimum value of the number of detection values included in the partial period.
  • the calculation unit 2040 may determine the number of feature constants ⁇ m ⁇ by a method described later, and may extract the number of partial periods from the domain of ⁇ g (t) ⁇ . In this case, for example, the calculating unit 2040 divides the domain of ⁇ g (t) ⁇ into a plurality of periods based on the magnitude of the gradient ⁇ g '(t) ⁇ , and selects the ⁇ m ⁇ The period is extracted as the above-described partial period.
  • ⁇ Feature Constant Calculation Method 2 In order to adopt the above-described calculation method 1, since it is necessary to take the log of y '(t), the value of y' (t) must always be positive. On the other hand, in the calculation method 2 described here, the value of y ′ (t) may be negative while using the index corresponding to the gradient g ′ (t) of g (t) described above. Can also be used.
  • the log of y ′ (t) is not taken, and y ′′ (t) / y ′ (t) is used as an index corresponding to the slope of g (t) g ′ (t).
  • Y '(t) may be negative.
  • g ′ (t) ⁇ can be said to be the direction of the vector (y ′ (t), ⁇ y ′′ (t)) ⁇ .
  • g ′ (t) is the direction of the velocity vector (vector representing time change) of the vector (y (t), ⁇ y ′ (t)) ⁇ .
  • the calculating unit 2040 calculates a vector ⁇ (y (t), ⁇ y '(t)) ⁇ for each time ⁇ t ⁇ , and further uses each calculated vector to calculate a speed vector ⁇ (y (t + 1)) for each time ⁇ t ⁇ .
  • -y (t), ⁇ y '(t + 1) -y' (t)) ⁇ is calculated.
  • the calculation unit 2040 extracts a plurality of partial periods in which the directions are substantially the same for the calculated plurality of velocity vectors from the measurement period of y (t).
  • y (t + 1) ⁇ is a detection value obtained after ⁇ y (t) ⁇ in the sensor.
  • the calculation unit 2040 then calculates the direction of the velocity vector of ⁇ (y (t), ⁇ y '(t)) ⁇ in each extracted partial period as a velocity constant corresponding to the partial period. For example, the calculating unit 2040 calculates a statistical value (such as an average value) of the direction of the speed vector in each partial period as a speed constant corresponding to the partial period. In addition, for example, the calculating unit 2040 performs, for each partial period, linear regression on the points ⁇ (y (t), ⁇ y ′ (t)) ⁇ included in the partial period, and calculates the slope of the regression line obtained for each partial period. May be a rate constant corresponding to the partial period.
  • the division of the measurement period of y (t) ⁇ ⁇ can be realized by clustering the velocity vectors in the directions. Note that the direction of the velocity vector ⁇ (y (t + 1) -y (t), ⁇ y '(t + 1) -y' (t)) ⁇ is ⁇ y '(t + 1) -y' (t) ⁇ / ⁇ y (t + 1) -y (t) ⁇ .
  • FIG. 8 is a diagram illustrating a velocity vector of (y (t), ⁇ y ′ (t)) ⁇ .
  • two partial periods are obtained from the time series data 14.
  • the calculation unit 2040 calculates the speed constant ⁇ 1 based on the direction of the speed vector in the partial period 1 and calculates the speed constant ⁇ 2 based on the direction of the speed vector in the partial period 2.
  • a direction different from any of the directions before and after the direction may appear for a short period of time. This is the same as the section 20 in FIG.
  • the direction of the velocity vector that appears only for such a short period may not be included in the feature constant.
  • a method of determining the minimum value of the length of the partial period or determining the number of characteristic constants can be used.
  • the contribution vector ⁇ can be calculated by performing parameter estimation on the contribution vector ⁇ using the time-series data 14 that is observation data.
  • An example of a prediction model in the case of using the rate constant ⁇ as the feature constant can be expressed by Expression (6).
  • an example of a prediction model when the time constant ⁇ is used as the feature constant can be expressed by Expression (7).
  • the calculation unit 2040 estimates the parameter ⁇ by the maximum likelihood estimation using the predicted value obtained from the prediction model and the observation value (that is, the time-series data 14) obtained from the sensor 10.
  • the maximum likelihood estimation for example, the least square method can be used.
  • the parameter ⁇ is determined according to the following objective function.
  • T represents the length of the time-series data 14 (the number of detected values).
  • y ⁇ (ti) represents a predicted value at time ti.
  • the vector ⁇ that minimizes the above objective function can be calculated using the following equation (11).
  • the vector Y (y (t0), y (t1), etc.
  • is a hyperparameter representing the weight given to the regularization term.
  • the parameter ⁇ can be determined according to the following equation (13).
  • is a hyperparameter and needs to be determined in advance.
  • the value of ⁇ is determined through test measurement and simulation.
  • it is preferable that the value of ⁇ is a small value so that the contribution value ⁇ does not vibrate.
  • FIG. 7 is a diagram exemplifying a feature amount obtained for a single molecule. From this figure, it can be seen that there is a trade-off between the blunting of the peak and the increase in the vibration depending on the value of ⁇ . Specifically, if ⁇ is too large, the vibration decreases, but the width of the peak increases. When the width of the peak is large, the result of measuring two molecules having similar rate constants appears as one large peak, so that it is difficult to distinguish these molecules. That is, the sensitivity is reduced. On the other hand, if ⁇ is too small, the peak width becomes small, but the vibration increases.
  • the purpose of the simulation is to evaluate the degree of occurrence of such peak dulling and vibration while changing ⁇ .
  • the contribution vectors ⁇ 1 and ⁇ 2 of two virtual single molecules having two different rate constants ⁇ 1 and ⁇ 2 are obtained. calculate. Then, the inner product of these two contribution vectors is calculated as follows.
  • This function ⁇ f ( ⁇ v) ⁇ attenuates while oscillating. Therefore, the width of the main lobe of this vibration can be quantified as “peak width” and the level of the side lobe can be quantified as “magnitude of vibration”. Therefore, ⁇ is determined by selecting a value of ⁇ such that the main lobe width is as narrow as possible and the level of the side lobe is as small as possible.
  • one of the advantages of suppressing the oscillation of the contribution value is that the feature value is robust against changes in the time constant and the speed constant, as described above. In other words, the feature value becomes robust against a temperature change.
  • the reason will be described.
  • the vibration of the contribution value is small, the distance between the contribution vectors before and after the translation is short. This means that when the time constant or the speed constant slightly changes, the feature amount also slightly changes. That is, it means that the robustness of the feature amount is high. Therefore, it can be said that suppressing the vibration of the feature amount improves the robustness of the feature amount.
  • the regularization in the least squares method is not limited to the above ⁇ L2 ⁇ regularization, and other regularizations such as L1 regularization may be introduced.
  • ⁇ ) and P ( ⁇ ) are defined by a multivariate normal distribution as follows, for example.
  • ⁇ , ⁇ ) is a multivariate normal distribution with mean ⁇ and covariance ⁇ .
  • ⁇ ⁇ 2 is a parameter representing the variance of the observation error.
  • is a covariance matrix of the prior distribution of ⁇ , and an arbitrary semi-definite matrix may be given in advance, or may be determined by a method described later.
  • ⁇ ) and P ( ⁇ ) may be determined by a Gaussian process (GP) as follows.
  • calculation unit 2040 may determine parameter ⁇ by Bayesian estimation using time-series data 14 that is an observed value. Specifically, the parameter ⁇ is determined by calculating the following conditional expected value.
  • Y] is a conditional expected value assuming that ⁇ and Y follow the probability distribution of Expression (18).
  • the feature vector ⁇ ⁇ maximizing the objective function (14) and the feature vector ⁇ obtained from the conditional expected value (19) can be both calculated by the following equation (20).
  • the hyperparameters set in advance include a) the form of the covariance function ⁇ ( ⁇ , ⁇ ′), b) the parameters of the covariance function, and c) the measurement error parameter ⁇ ⁇ 2. is there. Perform the following steps while changing these.
  • As an index for quantifying the magnitude and peak width of the vibration of the feature amount, for example, the width of the inlobe and the level of the side lobe of the function ⁇ f ( ⁇ v) ⁇ are used. In addition, a variance (square variance or absolute value variance) when the estimated ⁇ is regarded as a probability distribution may be used. These dispersion values become smaller as the vibration is smaller and the peak width is smaller. Note that an actual measurement (test measurement) may be performed instead of the simulation.
  • the calculating unit 2040 solves the combination optimization problem of minimizing the minimum value of the objective function of the least squares method with respect to the set of characteristic constants ⁇ instead of using the least squares method after calculating the characteristic constants.
  • a set of feature constants ⁇ and a set of contribution values ⁇ may be calculated. Specifically, the following objective function h ( ⁇ ) is minimized, and the parameter ⁇ ⁇ ⁇ ⁇ in the minimized h ( ⁇ ) is a calculation result of the calculation unit 2040.
  • an existing method can be used as a specific method for solving the above-described combination optimization problem.
  • the objective function becomes a monotonically decreasing dominant modular function, and thus can be accurately calculated by, for example, a greedy algorithm.
  • the calculation unit 2040 may determine the number of feature constants.
  • the determined number of characteristic constants can be used, for example, to determine the domain of g (t) and the number of partial periods extracted from the measurement period of y (t).
  • the number of feature constants can be determined using an information amount criterion such as AIC (Akaike's Information Criterion) or BIC (Bayesian Information Criterion).
  • AIC Kaike's Information Criterion
  • BIC Bayesian Information Criterion
  • AIC can be calculated as follows. The method of deriving the following AIC will be described later.
  • ⁇ ⁇ 2 is the variance of the observation error.
  • ⁇ Calculation unit 2040 determines integer ⁇ K ⁇ at which AIC ⁇ is minimized as the number of feature constants.
  • FIG. 9 is a diagram showing AIC in a graph. In FIG. 9, the upper part is a graph of the function ⁇ L ⁇ , and the lower part is a graph of ⁇ AIC ⁇ . The integer ⁇ K ⁇ that minimizes AIC is K *. Therefore, calculation unit 2040 sets the number of feature constants to ⁇ K * ⁇ .
  • the merit of improving the prediction accuracy obtained by increasing the number of parameters (here, the number of feature constants and contribution values) of the prediction model is obtained.
  • the number of feature constants can be appropriately determined in consideration of a balance with the disadvantage that the complexity of the model is increased by increasing the number of parameters.
  • a method for deriving the AIC of Expression (24) will be described.
  • the definition of AIC is as follows.
  • l is the maximum likelihood
  • p is the number of parameters.
  • the number of parameters p is 2K because the feature constant and the contribution value are each K.
  • the maximum likelihood is calculated as follows.
  • equation (24) By substituting the ⁇ l ⁇ and ⁇ p ⁇ into the definition of ⁇ AIC ⁇ , equation (24) can be obtained.
  • the output unit 2060 outputs information (hereinafter, output information) in which the set of feature constants ⁇ and the set of contribution values ⁇ obtained by the above-described method are associated with each other as a feature amount representing the feature of the gas (S106).
  • output information is text data representing a feature matrix F.
  • the output information may be information that graphically expresses the correspondence between the set of feature constants ⁇ and the set of contribution values ⁇ in a table, a graph, or the like.
  • FIG. 10 is a diagram exemplifying the characteristic amount of the target gas in a graph.
  • the horizontal axis represents the time constant ⁇
  • the vertical axis represents the contribution value ⁇ .
  • the pair of ( ⁇ i, ⁇ i) calculated by the calculation unit 2040 is plotted on the graph.
  • the output unit 2060 stores the output information in an arbitrary storage device.
  • the output unit 2060 causes the display device to display output information.
  • the output unit 2060 may transmit the output information to a device other than the information processing device 2000.
  • the information processing apparatus 2000 may calculate the association between the set of feature constants ⁇ ⁇ and the set of contribution values ⁇ for each of the plurality of time-series data 14 obtained for the same target gas.
  • the output unit 2060 may use a set of the plurality of associations as the feature amount of the target gas.
  • the information processing apparatus 2000 calculates a set of feature constants ⁇ u ⁇ and a set of contribution values ⁇ u ⁇ for the rising time-series data 14, and generates a feature matrix ⁇ Fu ⁇ in which these are correlated.
  • the information processing apparatus 2000 calculates a set of feature constants ⁇ d ⁇ and a set of contribution values ⁇ d ⁇ for the falling time-series data 14, and generates a feature matrix Fd ⁇ that associates these.
  • the information processing device 2000 outputs ⁇ Fu, ⁇ Fd ⁇ , which is a set of the generated feature matrices, as the feature amount of the target gas.
  • FIG. 11 is a diagram illustrating a case where a feature matrix is obtained from each of the rising time series data 14 and the falling time series data 14.
  • a feature matrix Fu is obtained from time-series data 14-1 having rising time-series data.
  • a feature matrix Fd is obtained from time-series data 14-2 which is time-series data at the falling edge. Therefore, the output unit 2060 outputs ⁇ Fu, ⁇ Fd ⁇ obtained by combining the obtained two feature matrices as the feature amount of the target gas.
  • the output unit 2060 outputs one matrix obtained by connecting the feature matrix obtained from the rising time-series data 14 and the feature matrix obtained from the falling time-series data 14 to the feature amount of the target gas. It may be.
  • the number of rows of the connected feature matrices is different from each other, the number of rows of the connected feature matrices is made to match each other by extending the matrix having the smaller number of rows by a method such as zero padding.
  • the plurality of feature matrices are not limited to those obtained from the rising time series data 14 and the falling time series data 14, respectively.
  • a plurality of time-series data 14 may be obtained by exposing each of the plurality of sensors 10 having different characteristics to the target gas.
  • the ease with which each molecule attaches to the sensor depends on the characteristics of the sensor.
  • the ease with which each molecule adheres to the functional film differs depending on the material of the functional film. The same applies to the ease of separation of each molecule. Therefore, by preparing the sensors 10 having the functional films of different materials and obtaining and analyzing the time-series data 14 from each of the plurality of sensors 10, the characteristics of the target gas can be grasped more accurately.
  • the information processing apparatus 2000 acquires the time-series data 14 from each of the plurality of sensors 10 having different characteristics, and generates information in which a set of feature constants and a set of contribution values are associated with each of the time-series data 14.
  • the output unit 2060 outputs a set of a plurality of pieces of information obtained as described above as a feature amount of the target gas.
  • FIG. 12 is a diagram illustrating a case where a plurality of feature matrices are obtained by obtaining the time-series data 14 from each of the plurality of sensors 10.
  • three sensors 10-1, 10-2, and 10-3 having different characteristics are prepared, and time-series data 14-1, time-series data 14-2, and time-series Data 14-3 is obtained.
  • the information processing device 2000 calculates the feature matrices F1, F2, and ⁇ F3 ⁇ from the plurality of time-series data 14, respectively. Then, the information processing apparatus 2000 outputs a set of these three feature matrices as the feature amount of the target gas. Note that, as described above, instead of outputting a set of a plurality of feature matrices, one feature matrix Fc obtained by connecting the plurality of feature matrices may be output.
  • the plurality of sensors 10 having different characteristics may be housed in one housing, or may be housed in different housings.
  • the sensor 10 is configured such that a plurality of functional films of different materials are housed in one sensor housing, and a detection value is obtained for each functional film.
  • the method described in FIG. 11 and the method described in FIG. 12 may be combined. That is, the information processing apparatus 2000 obtains the rising time-series data 14 and the falling time-series data 14 from each of the plurality of sensors 10, and calculates a feature matrix ⁇ F ⁇ for each of the obtained time-series data 14.
  • a set of a plurality of feature matrices described above or a single feature matrix obtained by connecting the sets may be used as the feature amount of the target gas.
  • the detection value of the sensor 10 may include a bias term that does not represent a change with time.
  • the time series data 14 is represented as follows.
  • a velocity constant ⁇ is used as the characteristic constant.
  • the bias is caused, for example, by the offset of the sensor 10 being shifted.
  • the bias is generated by the contribution of components commonly contained in the target gas and the purge gas (for example, the contribution of nitrogen and oxygen in the atmosphere).
  • the information processing device 2000 may have a function of removing an offset from the time-series data 14. This makes it possible to more accurately calculate the characteristic amount of the target gas.
  • a method of calculating the feature amount in consideration of the offset will be described.
  • the calculation unit 2040 calculates the contribution vector ⁇ ⁇ in consideration of the bias by generating the prediction model of the time-series data 14 represented by the above equation (27). That is, the calculation unit 2040 estimates the parameters ⁇ and ⁇ b ⁇ for the prediction model represented by Expression (27). Specifically, the calculation unit 2040 estimates ⁇ and ⁇ b ⁇ by optimizing the objective function (10), (12), or (16) not only for ⁇ but also for b. When a time constant is used as the feature constant, ⁇ k is replaced by 1 / ⁇ k in equation (27).
  • equation (14) is used as the objective function.
  • the calculation unit 2040 calculates ⁇ and b by the following optimization problem. The same applies when (8) or (10) is used as the objective function.
  • the output unit 2060 may output a bias ⁇ b ⁇ or ⁇ b0 ⁇ in addition to the feature matrix ⁇ F ⁇ . If the bias is due to the offset of the sensor offset, the value of b0 can be used to calibrate the sensor offset.
  • a time-series data acquisition unit that acquires the time-series data of the detection value, which is output from a sensor whose detection value changes in accordance with the attachment and detachment of molecules contained in the target gas, A plurality of feature constants that contribute to the time-series data, and a calculation unit that calculates a contribution value representing the magnitude of the contribution to the time-series data of each of the feature constants,
  • An output unit that outputs a combination of a plurality of characteristic constants and contribution values calculated for each characteristic constant as a characteristic amount of gas sensed by the sensor,
  • the information processing apparatus wherein the feature constant is a time constant or a rate constant relating to a magnitude of a temporal change in an amount of a molecule attached to the sensor.
  • the calculation unit Extracting a plurality of partial periods from the measurement period of the time-series data, For each of the partial periods, the characteristic constant is calculated based on a logarithm of a time change rate of the detection value in the partial period,
  • the partial period is a period in which the logarithm of the time change rate of the detection value included therein is substantially the same.
  • the calculation unit Using the time-series data, the detection value at each time and the time-series vector data with the time change rate of the detection value at that time as an element, A velocity vector is calculated for each of the calculated time-series vector data, Based on the direction of the velocity vector, from the measurement period of the time-series data, to extract a plurality of partial periods, For each of the partial periods, calculate the feature constant based on the direction of the velocity vector in the partial period, The partial period is a period in which directions of the velocity vectors included therein are substantially the same.
  • An information processing apparatus according to claim 1. 4.
  • the calculation unit calculates each contribution value by performing parameter estimation using the obtained time-series data for a prediction model of the detection value of the sensor using the contribution value of each of a plurality of feature constants as a parameter. 1. To 3. An information processing device according to any one of the above. 5. 3. The calculation unit calculates each of the contribution values by performing maximum likelihood estimation by a least squares method on the time series data obtained from the prediction model and the obtained time series data. An information processing apparatus according to claim 1. 6. 4. In the maximum likelihood estimation in the least square method, the objective function includes a regularization term. An information processing apparatus according to claim 1. 7. 3.
  • the calculation unit calculates each of the contribution values by MAP (Maximum a Posteriori) estimation or Bayes estimation using the prior distribution of each of the contribution values and the acquired time-series data.
  • the calculation unit for an objective function representing a square error between the time series data obtained from the prediction model and the acquired time series data, minimize the minimum value of the objective function with respect to a plurality of feature constants, 3. calculate the characteristic constant and a plurality of contribution values.
  • the prediction model includes a parameter representing a bias, 3. The calculation unit estimates parameters representing a contribution value and a bias for the prediction model.
  • the time series data acquisition unit acquires a plurality of time series data
  • the calculation unit calculates a set of a set of feature constants and a set of contribution values for each of the plurality of time-series data
  • the output unit outputs, as a feature amount of the target gas, information in which a plurality of sets of the calculated set of feature constants and the set of contribution values are collected.
  • An information processing device according to any one of the above. 12.
  • the plurality of time-series data includes both time-series data obtained when exposing the sensor to the target gas and time-series data obtained when removing the target gas from the sensor.
  • the plurality of time-series data includes time-series data obtained from each of the plurality of sensors having different characteristics. An information processing apparatus according to claim 1.
  • a control method executed by a computer A time-series data acquisition step of acquiring time-series data of the detection value, which is output from a sensor whose detection value changes in accordance with the attachment and detachment of molecules contained in the target gas, A plurality of feature constants that contribute to the time-series data, and a calculation step of calculating a contribution value representing the magnitude of the contribution of each of the feature constants to the time-series data, An output step of outputting a combination of a plurality of characteristic constants and contribution values calculated for each characteristic constant as a characteristic amount of the gas sensed by the sensor,
  • the control method wherein the characteristic constant is a time constant or a rate constant relating to a magnitude of a temporal change in an amount of a molecule attached to the sensor.
  • the characteristic constant is calculated based on a logarithm of a time change rate of the detection value in the partial period, 13.
  • the partial period is a period in which the logarithm of the temporal change rate of the detection value included therein is substantially the same. The control method described in 1. 16.
  • a velocity vector is calculated for each of the calculated time-series vector data, Based on the direction of the velocity vector, from the measurement period of the time-series data, to extract a plurality of partial periods, For each of the partial periods, calculate the feature constant based on the direction of the velocity vector in the partial period, 13.
  • the partial period is a period in which the directions of the velocity vectors included therein are substantially the same.
  • the contribution value is calculated by performing maximum likelihood estimation by a least squares method on the time series data obtained from the prediction model and the obtained time series data.
  • the maximum likelihood estimation in the least square method the objective function includes a regularization term. The control method described in 1. 20. 16.
  • the prior distribution is a multivariate normal distribution or a Gaussian process.
  • a minimum value of the objective function is minimized with respect to a plurality of feature constants. 17. calculating a characteristic constant and a plurality of contribution values of The control method described in 1. 23.
  • the prediction model includes a parameter representing a bias, 16.
  • the plurality of time-series data includes both time-series data obtained when exposing the sensor to the target gas and time-series data obtained when removing the target gas from the sensor.
  • the plurality of time-series data includes time-series data obtained from each of the plurality of sensors having different characteristics. The control method described in 1.

Abstract

情報処理装置(2000)は、センサ(10)によって出力された時系列データ(14)を取得し、複数の特徴定数θi と、各特徴定数θi の時系列データ(14)に対する寄与を表す寄与値ξi を算出する。そして情報処理装置(2000)は、特徴定数θi の集合Θと、寄与値ξi の集合Ξとを対応づけた情報を、対象ガスの特徴量として出力する。特徴定数θとしては、速度定数βや、速度定数の逆数である時定数τを採用することができる。

Description

情報処理装置、制御方法、及びプログラム
 本発明はガスの特徴の解析に関する。
 ガスをセンサで測定することにより、ガスに関する情報を得る技術が開発されている。特許文献1は、ナノメカニカルセンサで試料ガスを測定することで得られるシグナル(検出値の時系列データ)を利用して、試料ガスの種類を判別する技術を開示している。具体的には、センサの受容体に対する試料ガスの拡散時定数が、受容体の種類と試料ガスの種類の組み合わせによって決まるため、シグナルから得られる拡散時定数と、受容体の種類とに基づいて、試料ガスの種類を判別できることが開示されている。
特開2017-156254号公報
 特許文献1では、試料ガスに含まれている分子が1種類であることが前提となっており、複数種類の分子が混合している試料ガスを扱うことが想定されていない。本願発明は上記の課題に鑑みてなされたものであり、複数種類の分子が混合しているガスの特徴を抽出する技術を提供することである。
 本発明の情報処理装置は、1)対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得部と、2)時系列データに対して寄与する複数の特徴定数と、各特徴定数の時系列データに対する寄与の大きさを表す寄与値とを算出する算出部と、3)複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、センサによってセンシングされたガスの特徴量として出力する出力部と、を有する。特徴定数は、センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である。
 本発明の制御方法は、コンピュータによって実行される。当該制御方法は、1)対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得ステップと、2)時系列データに対して寄与する複数の特徴定数と、各特徴定数の時系列データに対する寄与の大きさを表す寄与値とを算出する算出ステップと、3)複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、センサによってセンシングされたガスの特徴量として出力する出力ステップと、を有する。特徴定数は、センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である。
 本発明のプログラムは、コンピュータに、本発明の制御方法が有する各ステップを実行させる。
 本発明によれば、複数種類の分子が混合しているガスの特徴を抽出する技術が提供される。
 上述した目的、およびその他の目的、特徴および利点は、以下に述べる好適な実施の形態、およびそれに付随する以下の図面によってさらに明らかになる。
実施形態1の情報処理装置の概要を例示する図である。 情報処理装置が取得するデータを得るためのセンサを例示する図である。 実施形態1の情報処理装置の機能構成を例示する図である。 情報処理装置を実現するための計算機を例示する図である。 実施形態1の情報処理装置によって実行される処理の流れを例示するフローチャートである。 センサから得られる複数の時系列データを例示する図である。 g(t) のグラフを例示する図である。 (y(t), y'(t)) の速度ベクトルを例示する図である。 AIC をグラフで表す図である。 対象ガスの特徴量をグラフで例示する図である。 立ち上がりの時系列データと立ち下がりの時系列データそれぞれから特徴行列を得るケースを例示する図である。 複数のセンサそれぞれから時系列データを得ることで、複数の特徴行列を得るケースを例示する図である。
 以下、本発明の実施の形態について、図面を用いて説明する。尚、すべての図面において、同様な構成要素には同様の符号を付し、適宜説明を省略する。また、特に説明する場合を除き、各ブロック図において、各ブロックは、ハードウエア単位の構成ではなく、機能単位の構成を表している。
[実施形態1]
<発明の概要と理論的背景>
 図1は、実施形態1の情報処理装置2000の概要を例示する図である。また、図2は、情報処理装置2000が取得するデータを得るためのセンサ10を例示する図である。センサ10は、分子が付着する受容体を有し、その受容体における分子の付着と離脱に応じて検出値が変化するセンサである。なお、センサ10によってセンシングされているガスを、対象ガスと呼ぶ。また、センサ10から出力される検出値の時系列データを、時系列データ14と呼ぶ。ここで、必要に応じ、時系列データ14を Y とも表記し、時刻 t の検出値を y(t) とも表記する。Y は、y(t) が列挙されたベクトルとなる。
 例えばセンサ10は、膜型表面応力(Membrane-type Surface Stress; MSS)センサである。MSS センサは、受容体として、分子が付着する官能膜を有しており、その官能膜に対する分子の付着と離脱によってその官能膜の支持部材に生じる応力が変化する。MSS センサは、この応力の変化に基づく検出値を出力する。なお、センサ10は、MSS センサには限定されず、受容体に対する分子の付着と離脱に応じて生じる、センサ10の部材の粘弾性や動力学特性(質量や慣性モーメントなど)に関連する物理量の変化に基づいて検出値を出力するものであればよく、カンチレバー式、膜型、光学式、ピエゾ、振動応答などの様々なタイプのセンサを採用することができる。
 ここで、説明のため、センサ10によるセンシングを以下のようにモデル化する。
(1)センサ10は、K 種類の分子を含む対象ガスに曝されている。
(2)対象ガスにおける各分子 k の濃度は一定のρkである。
(3)センサ10には、合計 N 個の分子が吸着可能である。
(4)時刻t においてセンサ10に付着している分子k の数は nk(t) 個である。
 センサ10に付着している分子 k の数 nk(t) の時間変化は、以下のように定式化できる。
Figure JPOXMLDOC01-appb-M000001
 式(1)の右辺の第1項と第2項はそれぞれ、単位時間当たりの分子 k の増加量(新たにセンサ10に付着する分子 k の数)と減少量(センサ10から離脱する分子 k の数)を表している。また、αk とβk はそれぞれ、分子 k がセンサ10に付着する速度を表す速度定数と、分子 k がセンサ10から離脱する速度を表す速度定数である。
 ここで、濃度ρkが一定であるため、上記式(1)から、時刻t における分子 k の数 nk(t) は、以下のように定式化できる。
Figure JPOXMLDOC01-appb-M000002
 また、時刻 t0(初期状態)でセンサ10に分子が付着していないと仮定すれば、nk(t) は以下のように表される。
Figure JPOXMLDOC01-appb-M000003
 センサ10の検出値は、対象ガスに含まれる分子によってセンサ10に働く応力によって定まる。そして、複数の分子によってセンサ10に働く応力は、個々の分子に働く応力の線形和で表すことができると考えられる。ただし、分子によって生じる応力は、分子の種類によって異なると考えられる。すなわち、センサ10の検出値に対する分子の寄与は、その分子の種類によって異なると言える。
 そこで、センサ10の検出値 y(t) は、以下のように定式化できる。
Figure JPOXMLDOC01-appb-M000004
 ここで、γk とξk はいずれも、センサ10の検出値に対する分子 k の寄与を表す。なお、「立ち上がり」と「立ち下がり」の意味については、後述する。
 ここで、対象ガスをセンシングしたセンサ10から得た時系列データ14を上述の式(4)のように分解できれば、対象ガスに含まれる分子の種類や、各種類の分子が対象ガスに含まれる割合を把握することができる。すなわち、式(4)に示す分解によって、対象ガスの特徴を表すデータ(すなわち、対象ガスの特徴量)が得られる。
 そこで情報処理装置2000は、センサ10によって出力された時系列データ14を取得し、時系列データ14を以下の式(5)に示すように分解する。
Figure JPOXMLDOC01-appb-M000005
 ここで、ξi は、センサ10の検出値に対する特徴定数θi の寄与を表す寄与値である。
 具体的には、まず情報処理装置2000は、時系列データ14を用いて、時系列データ14に対して寄与する複数の特徴定数の集合Θ={θ1,..., θm} と、各特徴定数θi の時系列データ14に対する寄与の大きさを表す寄与値ξi を算出する。なお、後述するように、特徴定数θi を算出してから寄与値ξi が算出されるケースと、特徴定数θi と共に寄与値ξi が算出されるケースとがある。
 さらに情報処理装置2000は、特徴定数の集合Θと寄与値の集合Ξとを対応づけた情報を、対象ガスの特徴を表す特徴量として出力する。特徴定数の集合Θと寄与値の集合Ξとの対応づけは、例えば、m 行2列の特徴行列 F で表される(m は特徴定数と寄与値それぞれの数)。例えばこの行列 F は、特徴定数の集合を表す特徴定数ベクトルΘ=(θ1,..., θm) を第1列に有し、なおかつ寄与値の集合を表す寄与ベクトルΞ=(ξ1,..., ξm) を第2列に有する。すなわち、F=(ΘT, ΞT) である。以降の説明では、特に断らない限り、対象ガスの特徴量は、特徴行列 F で表されるとする。ただし、対象ガスの特徴量は、必ずベクトルとして表現しなければならないわけではない。
 ここで、特徴定数θとしては、前述した速度定数βや、速度定数の逆数である時定数τを採用することができる。θとしてβとτを使う場合それぞれについて、式(5)は、以下のように表すことができる。
Figure JPOXMLDOC01-appb-M000006
 なお、図1では、特徴定数として速度定数βが利用されている。そのため、特徴定数の集合Θは、β1 からβm を示している。
<作用・効果>
 前述したように、センサ10の検出値に対する分子の寄与は、その分子の種類によって異なると考えられるため、上述した特徴定数の集合Θとそれに対応する寄与値の集合Ξは、対象ガスに含まれる分子の種類やその混合比率に応じて異なるものになると考えられる。よって、特徴定数の集合Θと寄与値の集合Ξとを対応づけた情報は、複数種類の分子が混合されているガスを互いに区別することができる情報、すなわちガスの特徴量として利用することができる。
 そこで本実施形態の情報処理装置2000は、対象ガスをセンサ10でセンシングすることで得られた時系列データ14に基づいて、特徴定数の集合Θと、各特徴定数の時系列データ14に対する寄与を表す寄与値の集合Ξとを算出し、算出した集合ΘとΞとを対応づけた情報を、対象ガスの特徴量として出力する。こうすることで、複数種類の分子が混合しているガスを識別することが可能な特徴量を、そのガスをセンサ10でセンシングした結果から自動的に生成することができる。
 なお、図1を参照した上述の説明は、情報処理装置2000の理解を容易にするための例示であり、情報処理装置2000の機能を限定するものではない。以下、本実施形態の情報処理装置2000についてさらに詳細に説明する。
<情報処理装置2000の機能構成の例>
 図3は、実施形態1の情報処理装置2000の機能構成を例示する図である。情報処理装置2000は、時系列データ取得部2020、算出部2040、及び出力部2060を有する。時系列データ取得部2020は、センサ10から時系列データ14を取得する。特徴定数生成部2030は、時系列データ14を用いて、複数の特徴定数θi と、各特徴定数θi の時系列データ14に対する寄与の大きさを表す寄与値ξi を算出する。出力部2060は、特徴定数の集合Θと寄与値の集合Ξとを対応づけた情報を、センサ10によってセンシングされたガスの特徴量として出力する。
<情報処理装置2000のハードウエア構成>
 情報処理装置2000の各機能構成部は、各機能構成部を実現するハードウエア(例:ハードワイヤードされた電子回路など)で実現されてもよいし、ハードウエアとソフトウエアとの組み合わせ(例:電子回路とそれを制御するプログラムの組み合わせなど)で実現されてもよい。以下、情報処理装置2000の各機能構成部がハードウエアとソフトウエアとの組み合わせで実現される場合について、さらに説明する。
 図4は、情報処理装置2000を実現するための計算機1000を例示する図である。計算機1000は任意の計算機である。例えば計算機1000は、Personal Computer(PC)やサーバマシンなどの据え置き型の計算機である。その他にも例えば、計算機1000は、スマートフォンやタブレット端末などの可搬型の計算機である。計算機1000は、情報処理装置2000を実現するために設計された専用の計算機であってもよいし、汎用の計算機であってもよい。
 計算機1000は、バス1020、プロセッサ1040、メモリ1060、ストレージデバイス1080、入出力インタフェース1100、及びネットワークインタフェース1120を有する。バス1020は、プロセッサ1040、メモリ1060、ストレージデバイス1080、入出力インタフェース1100、及びネットワークインタフェース1120が、相互にデータを送受信するためのデータ伝送路である。ただし、プロセッサ1040などを互いに接続する方法は、バス接続に限定されない。
 プロセッサ1040は、CPU(Central Processing Unit)、GPU(Graphics Processing Unit)、FPGA(Field-Programmable Gate Array)などの種々のプロセッサである。メモリ1060は、RAM(Random Access Memory)などを用いて実現される主記憶装置である。ストレージデバイス1080は、ハードディスク、SSD(Solid State Drive)、メモリカード、又は ROM(Read Only Memory)などを用いて実現される補助記憶装置である。
 入出力インタフェース1100は、計算機1000と入出力デバイスとを接続するためのインタフェースである。例えば入出力インタフェース1100には、キーボードなどの入力装置や、ディスプレイ装置などの出力装置が接続される。その他にも例えば、入出力インタフェース1100には、センサ10が接続される。ただし、センサ10は必ずしも計算機1000と直接接続されている必要はない。例えばセンサ10は、計算機1000と共有している記憶装置に時系列データ14を記憶させてもよい。
 ネットワークインタフェース1120は、計算機1000を通信網に接続するためのインタフェースである。この通信網は、例えば LAN(Local Area Network)や WAN(Wide Area Network)である。ネットワークインタフェース1120が通信網に接続する方法は、無線接続であってもよいし、有線接続であってもよい。
 ストレージデバイス1080は、情報処理装置2000の各機能構成部を実現するプログラムモジュールを記憶している。プロセッサ1040は、これら各プログラムモジュールをメモリ1060に読み出して実行することで、各プログラムモジュールに対応する機能を実現する。
<処理の流れ>
 図5は、実施形態1の情報処理装置2000によって実行される処理の流れを例示するフローチャートである。時系列データ取得部2020は時系列データ14を取得する(S102)。算出部2040は、時系列データ14を用いて、複数の特徴定数θi と、各特徴定数θi に対応する寄与値ξiとを算出する(S104)。出力部2060は、特徴定数の集合と寄与値の集合とを対応づけた情報を、対象ガスの特徴量として出力する(S106)。
 情報処理装置2000が図5に示す一連の処理を実行するタイミングは様々である。例えば、情報処理装置2000は、時系列データ14を指定する入力操作を受け付け、指定された時系列データ14について一連の処理を実行する。その他にも例えば、情報処理装置2000は、時系列データ14を受信できるように待機しておき、時系列データ14を受信したこと(すなわち、S102を実行したこと)に応じて、S104以降の処理を実行する。
<時系列データ14の取得:S102>
 時系列データ取得部2020は時系列データ14を取得する(S102)。時系列データ取得部2020が時系列データ14を取得する方法は任意である。例えば情報処理装置2000は、時系列データ14が記憶されている記憶装置にアクセスすることで、時系列データ14を取得する。時系列データ14が記憶されている記憶装置は、センサ10の内部に設けられていてもよいし、センサ10の外部に設けられていてもよい。その他にも例えば、時系列データ取得部2020は、センサ10から出力される検出値を順次受信することで、時系列データ14を得てもよい。
 時系列データ14は、センサ10が出力した検出値を、センサ10から出力された時刻が早い順に並べた時系列のデータである。ただし、時系列データ14は、センサ10から得られた検出値の時系列データに対して、所定の前処理が加えられたものであってもよい。また、前処理が行われた時系列データ14を取得する代わりに、時系列データ取得部2020が時系列データ14に対して前処理を行ってもよい。前処理としては、例えば、時系列のデータからノイズ成分を除去するフィルタリングなどを採用することができる。
 ここで、時系列データ14は、センサ10を対象ガスに曝すことで得られる。ただし、センサを用いてガスに関する測定を行う場合、センサを測定対象のガスに曝す操作と、センサから測定対象のガスを取り除く操作を繰り返すことで、センサから解析対象の時系列データを複数得ることがある。
 図6は、センサから得られる複数の時系列データを例示する図である。図6では、立ち上がりの時系列データと立ち下がりの時系列データとを区別しやすいように、立ち上がりの時系列データが実線で表されており、立ち下がりの時系列データが点線で表されている。図6において、期間P1の時系列データ14-1と期間P3の時系列データ14-3は、センサを測定対象のガスに曝す操作によって得られる。このように、センサを測定対象のガスに曝す操作を行うと、センサの検出値が増加する。そこで、センサを測定対象のガスに曝す操作によって得られる時系列データを、「立ち上がり」の時系列データとも呼ぶ。式(4)における「立ち上がりの場合」は、「時系列データ14が立ち上がりの時系列データである場合」を意味する。以降の式においても同様である。
 一方、期間P2の時系列データ14-2と期間P4の時系列データ14-4は、センサから測定対象のガスを取り除く操作によって得られる。なお、センサから測定対象のガスを取り除く操作は、例えば、センサをパージガスと呼ばれるガスに曝すことで実現される。このように、センサから測定対象のガスを取り除く操作を行うと、センサの検出値が減少する。そこで、センサから測定対象のガスを取り除く操作によって得られる時系列データを、「立ち下がり」の時系列データとも呼ぶ。式(4)における「立ち下がりの場合」は、「時系列データ14が立ち下がりの時系列データである場合」を意味する。以降の式においても同様である。
 情報処理装置2000では、センサ10を対象ガスに曝す操作とセンサ10から対象ガスを取り除く操作のそれぞれで得られる時系列データ14とが区別され、それぞれ異なる時系列データ14として扱われる。例えば図6の例では、期間P1からP4という4つの期間それぞれで得られる時系列データが、それぞれ異なる時系列データ14として扱われる。そのため、センサ10を対象ガスに曝す操作とセンサ10から対象ガスを取り除く操作とを繰り返すことで一連の時系列データを得た場合には、その一連の時系列データを複数の時系列データ14に分割する必要がある。
 センサ10から得られた一連の時系列データを分割して複数の時系列データ14を得る方法には、様々な方法を採用することができる。例えば、センサ10から得られた一連の時系列データを人手で分割することで、複数の時系列データ14を得る。その他にも例えば、情報処理装置2000が一連の時系列データを取得し、その時系列データを分割して複数の時系列データ14を得てもよい。
 なお、情報処理装置2000によって時系列データを分割する方法には、様々な方法を採用できる。例えば、以下のような方法がある。
<<(a)1階微分を用いる方法>>
 時系列データ14において、分割すべき部分では、センサ値の微分が不連続となり、直後に絶対値が最大となる。そこで、1階微分の絶対値が大きくなる点を用いて、時系列データ14を分割することができる。
<<(b)2階微分を用いる方法>>
 同様に、分割すべき点では微分が不連続となるため、2階微分は無限大に発散する。したがって、2階微分の絶対値が大きくなる点を用いて、時系列データ14を分割することができる。
<<(c)センサから得られるメタデータを利用する方法>>
 センサの種類によっては、検出値以外のメタデータが提供される。例えば MSS のモジュールでは、測定対象のガス(サンプル)とパージガスの吸引にそれぞれ異なるポンプ(サンプルポンプとパージポンプ)が用意されており、これらを交互にオン/オフすることで、立ち上がりの測定と立ち下がりの測定が行われている。そして、記録される検出値には、ポンプの動作シーケンス(どちらのポンプを利用して得られた検出値であるかを表す情報や、流量のフィードバック制御に用いる流量計測値など)が、時系列情報として付加される。そこで例えば、情報処理装置2000は、時系列データ14と共に得られるポンプの動作シーケンスを利用することで、時系列データ14を分割することができる。
<<上記方法の組み合わせ>>
 (c)の方法では、ポンプが動作してからセンサにガスが届くまでの遅延を加味した修正を加えることが好ましい。そこで例えば、情報処理装置2000は、(c)の方法で時系列データ14を仮に複数の区間に分割した後、各区間の中で1階微分の絶対値が最大となる時点を特定し、特定した各時点で時系列データ14を分割する。
 なお、情報処理装置2000は、センサ10を対象ガスに曝す操作で得られる時系列データ14と、センサ10から対象ガスを取り除く操作で得られる時系列データ14のうち、いずれか一方のみを利用するように構成されてもよい。
<特徴定数と寄与値の算出:S104>
 算出部2040は、時系列データ14を用いて、複数の特徴定数θi と、特徴定数θi に対応する寄与値ξi とを算出する(S104)。この処理は、時系列データ14を、式(5)に示したξi*f(θi) の和に分解することに相当する。
 算出部2040が時系列データ14から特徴定数とその寄与値を算出する方法には、様々な方法を利用することができる。以下、まずは特徴定数の算出方法について説明する。なお、特徴定数を算出する過程でその寄与値が算出されるケースもある。このケースについては後述する。
<<特徴定数の算出方法1>>
 例えば算出部2040は、時系列データ14の微分 y'(t) に log をとることで得られる log[y'(t)] の傾きに基づいて、特徴定数の集合を算出する。log[y'(t)] は、以下のようになる。
Figure JPOXMLDOC01-appb-M000007
 なお、log をとるため、y'(t) が常に正であることを前提とする。y'(t) が負になるケースでは、後述する他の方法を利用して、特徴定数を算出する。
 図7は、g(t) のグラフを例示する図である。式(8)より、g(t)=log[y'(t)] は、いわゆる log-sum-exp と呼ばれる形式の関数であり、この関数は、積算されている exp のうちの最大値に近い値をとることが知られている。すなわち、式(8)では、g(t)の値が、その時刻 t における ci-βit のうちで最大のものに近い値となる。
 ここで、近似的に g(t)=ci-βit とおくと、この関数の傾き g'(t) は速度定数 βi となる。よって、時刻 t における g(t) の傾き g'(t) は、その時刻において ci-βit が最大となる i に対応するβi と略同一になる。また、exp(ci-βit) は単調減少関数であることから、図7に示すように、g(t) の傾きは段階的に切り替わると言える。
 そこで算出部2040は、各時刻 t における g(t) の傾き g'(t) を算出し、g(t) の定義域(時系列データ14の測定期間)から、g'(t) が略同一である複数の期間(以下、部分期間)を複数抽出する。そして算出部2040は、抽出した各部分期間における g'(t) を、速度定数βi として算出する。例えば算出部2040は、各部分期間における g'(t) の統計値(平均値など)を、その部分期間に対応する速度定数として算出する。その他にも例えば、算出部2040は、部分期間ごとに、その部分期間に含まれる (t, g(t)) について線形回帰を行い、各部分期間について得られる回帰直線の傾きを、速度定数としてもよい。なお、g'(t) が略同一である期間は、例えば、g'(t) の集合をその値の大きさでクラスタリングすることで特定することができる。
 例えば図7において、g(t) の定義域から、3つの部分期間が抽出されている。算出部2040は、これら3つの部分期間それぞれに対応する速度定数として、β1、β2、及びβ3 を算出する。このようにして、特徴定数の集合Θ={β1, β2, β3} が算出される。
 ここで、図7に示す区間20のように、g(t) の傾きが切り替わる部分に、その前後の傾きのいずれとも異なる傾きが、短い期間現れることがある。このような短い期間だけ現れる傾きについては、特徴定数に含めないようにしてもよい。
 例えば、部分期間の長さの最小値を予め定めておく。そして、算出部2040は、傾き g'(t) の大きさに基づいて g(t) の定義域を複数の期間に分割し、最小値以上の長さを持つ期間のみを、前述した部分期間として抽出する(すなわち、最小値未満の長さの期間に対応する g'(t) は特徴定数の集合に含めない)。こうすることで、短い期間だけ現れる g'(t) については、特徴定数に含めないようにすることができる。なお、部分期間の長さの最小値は、部分期間に含まれる検出値の数の最小値として表すこともできる。
 その他にも例えば、算出部2040は、後述する方法で特徴定数の個数 m を決定し、その個数の部分期間を g(t) の定義域から抽出してもよい。この場合、例えば算出部2040は、傾き g'(t) の大きさに基づいて g(t) の定義域を複数の期間に分割し、その中から期間の長さについて上位にある m 個の期間を、前述した部分期間として抽出する。
<<特徴定数の算出方法2>>
 前述した算出方法1を採用するためには、y'(t) の log をとる必要があることから、y'(t) の値が常に正である必要がある。これに対し、ここで説明する算出方法2では、前述した g(t) の傾き g'(t) に相当する指標を利用しつつ、y'(t) の値が負になることがあっても利用することができる。
 まず、前述した g(t)=log[y'(t)] の傾きである g'(t) は、以下のように表すことができる。
Figure JPOXMLDOC01-appb-M000008
 そこで算出方法2では、y'(t) の log をとらず、y''(t)/y'(t) を g(t) の傾き g'(t) に相当する指標として利用することで、y'(t) が負になることがあるケースも扱えるようにする。
 ここで、g'(t) は、(y'(t), y''(t)) というベクトルの向きということができる。さらに言い換えれば、g'(t) は、(y(t), y'(t)) というベクトルの速度ベクトル(時間変化を表すベクトル)の向きであるといえる。
 そこで算出部2040は、各時刻 t についてベクトル (y(t), y'(t)) を算出し、さらに算出した各ベクトルを利用して、各時刻 t について速度ベクトル (y(t+1)-y(t), y'(t+1)-y'(t)) を算出する。そして算出部2040は、算出された複数の速度ベクトルについて、その向きが略同一である複数の部分期間を、y(t) の測定期間から抽出する。なお、y(t+1) は、センサにおいて y(t) の次に得られる検出値である。
 そして、算出部2040は、抽出した各部分期間における (y(t), y'(t)) の速度ベクトルの向きを、その部分期間に対応する速度定数として算出する。例えば算出部2040は、各部分期間における上記速度ベクトルの向きの統計値(平均値など)を、その部分期間に対応する速度定数として算出する。その他にも例えば、算出部2040は、部分期間ごとに、その部分期間に含まれる点 (y(t), y'(t)) について線形回帰を行い、各部分期間について得られる回帰直線の傾きを、その部分期間に対応する速度定数としてもよい。なお、y(t) の測定期間の分割は、上記速度ベクトルをその向きでクラスタリングすることで実現することができる。なお、速度ベクトル (y(t+1)-y(t), y'(t+1)-y'(t)) の向きは、{y'(t+1)-y'(t)}/{y(t+1)-y(t)} で表される。
 図8は、(y(t), y'(t)) の速度ベクトルを例示する図である。図8では、時系列データ14から2つの部分期間が得られている。算出部2040は、部分期間1における速度ベクトルの向きに基づいて速度定数β1 を算出し、部分期間2における速度ベクトルの向きに基づいて速度定数β2 を算出している。
 ここで、図8に示す区間30のように、速度ベクトルの向きが変化する部分に、その前後の向きのいずれとも異なる向きが、短い期間現れることがある。これは、図7の区間20と同様である。このような短い期間だけ現れる速度ベクトルの向きについては、特徴定数に含めないようにしてもよい。その具体的な方法には、方法1で説明したように、部分期間の長さの最小値を定めておいたり、特徴定数の個数を決定したりする方法を利用することができる。
<<寄与値の算出方法>>
 ここでは、特徴定数を算出した後に、その特徴定数に対応する寄与値を算出する方法を説明する。算出部2040は、算出された特徴定数の集合Θ={θ1,..., θm} に対応する寄与値ξiの集合(すなわち、寄与ベクトル)Ξ={ξ1,..., ξm} をパラメータとして、センサ10の検出値を予測する予測モデルを生成する。この予測モデル生成する際、観測データである時系列データ14を利用して寄与ベクトルΞについてパラメータ推定を行うことにより、寄与ベクトルΞを算出することができる。特徴定数として速度定数βを使う場合の予測モデルの一例は、式(6)で表すことができる。また、特徴定数として時定数τを使う場合の予測モデルの一例は、式(7)で表すことができる。
 予測モデルのパラメータ推定には、種々の方法を利用することができる。以下、その方法についていくつか例示する。なお、以下の説明では、速度定数βを特徴定数として利用するケースに説明している。時定数τを特徴定数とする場合におけるパラメータ推定の方法は、以下の説明における速度定数βを1/τと読み替えることで実現できる。
<<パラメータ推定の方法1>>
 例えば算出部2040は、予測モデルから得られる予測値と、センサ10から得られた観測値(すなわち、時系列データ14)とを用いた最尤推定により、パラメータΞを推定する。最尤推定には、例えば最小二乗法を用いることができる。この場合、具体的には、以下の目的関数に従ってパラメータΞを決定する。
Figure JPOXMLDOC01-appb-M000009
 ここで、T は時系列データ14の長さ(検出値の数)を表す。また、y^(ti) は、時刻 ti の予測値を表す。
 上述の目的関数を最小化するベクトルΞは、以下の式(11)を用いて算出することができる。
Figure JPOXMLDOC01-appb-M000010
 ここで、ベクトル Y=(y(t0), y(t1),...) である。
 そこで、算出部2040は、時系列データ Yと特徴定数の集合Θ={β1, β2,...}を上記式(11)に適用することで、パラメータΞを算出する。
<<パラメータ推定の方法2>>
 上述した最小二乗法において、正則化項を導入して正則化を行ってもよい。例えば以下の式(12)は、L2 正則化を行う例を示している。 
Figure JPOXMLDOC01-appb-M000011
 ここで、λは正則化項に与える重みを表すハイパーパラメータである。
 この場合、以下の式(13)に従ってパラメータΞを決定することができる。
Figure JPOXMLDOC01-appb-M000012
 このような正則化項を導入すると、正則化項を導入しない場合と比較し、行列計算において測定誤差が増幅されてしまうことを抑えることができるため、各寄与値ξi をより正確に算出することができる。また、誤差の増幅を抑えることで、寄与値ξが数値的に安定するため、混合比に対する特徴量のロバスト性が向上する。
 なお、前述したように、λはハイパーパラメータであり、予め定めておく必要がある。例えば、テスト測定やシミュレーションを通じて、λの値を決定する。ここで、λの値は、寄与値ξが振動しない程度に小さい値とすることが好ましい。
 ここでは、λの値を決定するためのシミュレーションについて説明する。シミュレーションでは、仮想的に「寄与が1の単一の分子」を測定した場合(たとえば、立下りの場合、単一の分子の速度定数をβ0とすれば、y(t)=exp{-β0*t} となる)を考え、この場合の式(13)による特徴量推定値の結果を観察する。仮想的に、理想的な観測(無限小の測定間隔で無限長時間の測定ができ、観測誤差がゼロ)が可能であるとした場合、仮想的な単一分子のシミュレーションでは、以下のようにβ0だけに尖ったピークがある特徴量が得られ、もとの速度定数β=β0と寄与ξ=1 が完全に再現されるはずである。
Figure JPOXMLDOC01-appb-M000013
 しかし、実際には理想的な観測は不可能であるため、寄与値のピークが鈍化してしまったり、寄与値が振動してしまったりする。図7は、単一分子について得られる特徴量を例示する図である。この図から、λの値による、ピークの鈍化と振動の増加のトレードオフが分かる。具体的には、λが大きすぎると、振動は減るものの、ピークの幅が大きくなってしまう。ピークの幅が大きくなると、速度定数が近い2つの分子を測定した結果が、1つの大きなピークに見えてしまうため、これらの分子を見分けることが難しくなる。すなわち、感度が下がることになる。一方で、λが小さすぎると、ピークの幅は小さくなるものの、振動が増えてしまう。振動が増えると、後述するように、特徴量のロバスト性が低下してしまう。よって、振動が発生しない(ロバスト性を損なわない)程度にλを小さくするように決定することで、ピークを鋭く(感度を向上)することが好ましいと言える。
 シミュレーションの目的は、λを変化させながら、このようなピークの鈍化や振動の発生度合を評価することである。「振動の大きさ」と「ピーク幅」を定量的に測るためには、例えばシミュレーションにより、2つの異なる速度定数β1、β2をそれぞれ持つ2つの仮想的な単一分子の寄与ベクトルΞ1、Ξ2を計算する。そして、これら2つの寄与ベクトルの内積を以下のように計算する。
Figure JPOXMLDOC01-appb-M000014
 この関数 f(Δv) は、振動しながら減衰する。そこで、この振動のメインローブの幅を「ピーク幅」、サイドローブのレベルを「振動の大きさ」として定量化することができる。そこで、メインローブ幅がなるべく細く、なおかつサイドローブのレベルがなるべく小さくなるようなλの値を選択することにより、λを決定する。
 ここで、寄与値の振動を抑えることの利点の1つは、前述したように、特徴量が時定数や速度定数の変化に対してロバストになることである。言い換えれば、特徴量が温度変化に対してロバストになる。以下、その理由を説明する。
 温度変化によって時定数や速度定数に変化が生じると、図7や後述する図8に示す寄与値は、X軸方向に平行移動することになる。寄与値が大きく振動していると、寄与値がX軸方向に少し平行移動しただけでも、平行移動前後の寄与ベクトルの距離が大きくなってしまう。すなわち、時定数や速度定数が少し変化しただけで特徴量が大きく変化してしまい、時定数変化や速度定数変化に対する特徴量のロバスト性が低くなる。
 これに対し、寄与値の振動が少なければ、平行移動前後の寄与ベクトルの距離は短くなる。これは、時定数や速度定数が少し変化した場合に、特徴量も少しだけ変化することを意味する。すなわち、特徴量のロバスト性が高いことを意味する。よって、特徴量の振動を抑えることで、特徴量のロバスト性が向上すると言える。
 なお、最小二乗法における正則化は、前述した L2 正則化には限定されず、L1 正則化などの他の正則化を導入してもよい。
<<パラメータ推定の方法3>>
 この方法では、パラメータΞについて事前分布 P(Ξ) を設定しておく。そして算出部2040は、観測値である時系列データ14を用いた MAP (Maximum a Posteriori)推定により、パラメータΞを決定する。具体的には、以下の目的関数を最大化するパラメータΞを採用する。
Figure JPOXMLDOC01-appb-M000015
 P(Y|Ξ)とP(Ξ)は、例えば次のように多変量正規分布で定義される。
Figure JPOXMLDOC01-appb-M000016
 ここで、N(・|μ,Σ)は、平均μ、共分散Σの多変量正規分布である。また、ベクトル y^=(y^(t1), y^(t2),...)=ΦΞである。σ^2 は観測誤差の分散を表すパラメータである。
 ΛはΞの事前分布の共分散行列であり、予め任意の半正定値行列を与えても良いし、後述の方法などにより決定しても良い。
 また、P(Y|Ξ)、P(Ξ)は、次のように、ガウス過程(Gaussian Process; GP)によって定めても良い。
Figure JPOXMLDOC01-appb-M000017
 ここで、GP(ξ(β)|μ(β),Λ(β,β’))は、平均値関数がμ(β)、共分散関数(カーネル関数)がΛ(β,β’)のガウス過程である。また、ガウス過程は連続関数を生成する確率過程であるため、ここでは、ξ(β)は、寄与率をβ(もしくはτ)に関して表した連続関数であり、ベクトルΞは、関数ξ(β)の「β=β1,β2,・・・」における値を配列したベクトルΞ=(ξ(β1),ξ(β2),...)である。この場合、式(17) は式(18)の特別な場合とみなすことができ、式(17) における共分散行列Λの(i,j)成分は、式(18)の共分散関数Λ(β,β’)の(β,β’)=(β1, β2)における値である。すなわち、式(17)における行列Λは、所謂ガウス過程におけるグラム行列である。
 また、算出部2040は、観測値である時系列データ14を用いたベイズ推定により、パラメータΞを決定しても良い。具体的には、以下の条件付き期待値を計算することにより、パラメータΞを決定する。
Figure JPOXMLDOC01-appb-M000018
 ここで、E[Ξ|Y] は、ΞおよびYが式(18)の確率分布に従っていると仮定した場合の、条件付き期待値である。
 上記目的関数(14)を最大化する特徴ベクトルΞ、および、上記条件付き期待値(19)によって求まる特徴ベクトルΞは、いずれも以下の式(20)によって算出することができる。
Figure JPOXMLDOC01-appb-M000019
<<<ハイパーパラメータの定め方>>>
 ガウス過程を利用する場合、事前に設定しておくハイパーパラメータとして、a)共分散関数 Λ(β,β’)の形、b)共分散関数のパラメータ、及びc)測定誤差パラメータ σ^2 がある。これらを変えながら、次の手順を実行する。
1.仮想的な速度定数β0の単一分子の測定値をシミュレーションする。
2.シミュレーションした測定値から、寄与値を推定する。
3.推定した寄与値の振動の大きさ・ピーク幅を定量化する。
4.上述の a ~ c のハイパーパラメータを変更しながら、1~3を繰り返す。
5.グリッドサーチや最急降下法により、振動が小さく、ピーク幅が狭くなるように a ~ c のハイパーパラメータを決定する。
 なお、特徴量の振動の大きさとピーク幅を定量化する指標には、例えば、上述の関数 f(Δv) のインローブの幅とサイドローブのレベルを用いる。また、そのほかにも、推定したΞを確率分布として見做した際の分散値(二乗分散や絶対値分散)を用いてもよい。これらの分散値は、振動が小さく、ピーク幅が狭くなるほど、小さい値となる。なお、シミュレーションの代わりに、実際の測定(テスト測定)を実施してもよい。
<特徴定数と共に寄与値が算出されるケース>
 前述したように、特徴定数に対応する寄与値の算出方法の1つとして、最小二乗法を利用する方法がある。算出部2040は、特徴定数を算出した後に最小二乗法を利用する代わりに、前述した最小二乗法の目的関数の最小値を特徴定数の集合Θに関して最小化するという組み合わせ最適化問題を解くことにより、特徴定数の集合Θと寄与値の集合Ξの双方を算出してもよい。具体的には、以下の目的関数 h(Θ) を最小化するΘと、最小化された h(Θ) におけるパラメータΞが、算出部2040の算出結果となる。
Figure JPOXMLDOC01-appb-M000020
 ここで、上記の組み合わせ最適化問題を解く具体的な方法には、既存の方法を利用することができる。特に、Ξの要素がすべて正の数であるとき、上記の目的関数は単調減少な優モジュラ関数になるため、たとえば、貪欲アルゴリズムによって精度よく計算することができる。
<特徴定数の個数を決定する方法>
 算出部2040は、特徴定数の個数を決定してもよい。決定した特徴定数の個数は、例えば前述した、g(t) の定義域や y(t) の測定期間から抽出する部分期間の数を定めるために利用することができる。
 例えば特徴定数の個数は、AIC(Akaike's Information Criterion)や BIC(Bayesian Information Criterion)などの情報量基準を利用して決定することができる。情報量規準の算出には、前述した式(22)の関数 h(Θ) を利用する。例えば AIC は、以下の様に算出することができる。なお、下記 AIC の導出方法については後述する。
Figure JPOXMLDOC01-appb-M000021
 ここで、σ^2 は観測誤差の分散である。
 算出部2040は、AIC が最小となる整数 K を、特徴定数の個数として決定する。図9は、AIC をグラフで表す図である。図9において、上段は関数 L のグラフであり、下段は AIC のグラフである。AIC が最小になる整数 K は、K* である。そのため、算出部2040は、特徴定数の個数を K* とする。
 このように情報量規準を利用して特徴定数の個数を決定することにより、予測モデルのパラメータ(ここでは特徴定数や寄与値の数)の個数を増やすことによって得られる予測精度の向上というメリットと、パラメータの個数を増やすことによってモデルの複雑さが増すというデメリットとのバランスを考慮して、特徴定数の個数を適切に決定することができる。
 式(24)の AIC を導出する方法について説明する。まず、AIC の定義は以下の通りである。
Figure JPOXMLDOC01-appb-M000022
ただし、l は最大尤度であり、p はパラメータ数である。
 パラメータ数 p は、特徴定数と寄与値が K 個ずつであるため、2K である。また、最大尤度は以下のように算出される。
Figure JPOXMLDOC01-appb-M000023
 上記 l と p を AIC の定義に代入することで、式(24)を得ることができる。
<特徴量の出力:S106>
 出力部2060は、前述した方法で得られた特徴定数の集合Θと寄与値の集合Ξとを対応づけた情報(以下、出力情報)を、ガスの特徴を表す特徴量として出力する(S106)。例えば出力情報は、特徴行列 F を表すテキストデータである。その他にも例えば、出力情報は、特徴定数の集合Θと寄与値の集合Ξとの対応付けを、表やグラフなどでグラフィカルに表現した情報であってもよい。
 図10は、対象ガスの特徴量をグラフで例示する図である。図10のグラフは、横軸に時定数τを示し、縦軸に寄与値ξを示している。そして、算出部2040によって算出された(τi, ξi)のペアがグラフ上にプロットされている。このように対象ガスの特徴量をグラフィカルな情報で表すことにより、人がガスの特徴を直感的に理解しやすくなる。
 出力情報を出力する具体的な方法は様々である。例えば出力部2060は、出力情報を任意の記憶装置に記憶させる。その他にも例えば、出力部2060は、出力情報をディスプレイ装置に表示させる。その他に例えば、出力部2060は、情報処理装置2000以外の装置に出力情報を送信してもよい。
<特徴定数の集合と寄与値の集合との対応付けを複数算出するケース>
 情報処理装置2000は、同一の対象ガスについて得られた複数の時系列データ14それぞれについて、特徴定数の集合Θと寄与値の集合Ξとの対応付けを算出してもよい。そしてこの場合、出力部2060は、これら複数の対応付けの集合を、対象ガスの特徴量としてもよい。
 例えば情報処理装置2000は、立ち上がりの時系列データ14について、特徴定数の集合Θu と寄与値の集合Ξu とを算出し、これらを対応づけた特徴行列 Fu を生成する。また、情報処理装置2000は、立ち下がりの時系列データ14について、特徴定数の集合Θd と寄与値の集合Ξd とを算出し、これらを対応づけた特徴行列Fd を生成する。そして情報処理装置2000は、生成した特徴行列の組である {Fu, Fd} を対象ガスの特徴量として出力する。
 図11は、立ち上がりの時系列データ14と立ち下がりの時系列データ14それぞれから特徴行列を得るケースを例示する図である。図11において、立ち上がりの時系列データある時系列データ14-1から、特徴行列Fu が得られている。また、立ち下がりの時系列データである時系列データ14-2から、特徴行列Fd が得られている。そこで出力部2060は、得られた2つの特徴行列を組み合わせた {Fu, Fd} を、対象ガスの特徴量として出力する。
 なお、出力部2060は、立ち上がりの時系列データ14から得られる特徴行列と、立ち下がりの時系列データ14から得られる特徴行列とを連結することで得られる1つの行列を、対象ガスの特徴量としてもよい。例えばこの場合、出力部2060は、Fu=(ΘuT, ΞuT) とFd=(ΘdT, ΞdT) を連結したFc=(ΘuT, ΞuT, ΘdT, ΞdT) を、対象ガスの特徴量として出力する。なお、連結する特徴行列の行数が互いに異なる場合、行数が少ない方の行列をゼロパディング等の方法で拡張することで、連結する特徴行列の行数を互いに一致させる。
 複数の特徴行列は、立ち上がりの時系列データ14と立ち下がりの時系列データ14それぞれから得られるものに限定されない。例えば、特性の異なる複数のセンサ10それぞれを対象ガスに曝すことで、複数の時系列データ14を得てもよい。分子をセンサに付着させる場合、センサに対する各分子の付着しやすさは、センサの特性によって異なる。例えば官能膜に分子が付着するタイプのセンサを利用する場合、官能膜の材質によって、その官能膜に対する各分子の付着しやすさが異なる。各分子の離脱しやすさについても同様である。そのため、それぞれ異なる材質の官能膜を持つセンサ10を用意し、これら複数のセンサ10それぞれから時系列データ14を得て解析することで、対象ガスの特徴をより正確に把握することができる。
 そこで情報処理装置2000は、特性の異なる複数のセンサ10それぞれから時系列データ14を取得し、各時系列データ14について、特徴定数の集合と寄与値の集合とを対応づけた情報を生成する。出力部2060は、このようにして得られた複数の情報の組を、対象ガスの特徴量として出力する。
 図12は、複数のセンサ10それぞれから時系列データ14を得ることで、複数の特徴行列を得るケースを例示する図である。この例では、それぞれ特性が異なる3つのセンサ10-1、センサ10-2、及びセンサ10-3が用意されており、それぞれから時系列データ14-1、時系列データ14-2、及び時系列データ14-3が得られている。情報処理装置2000は、これら複数の時系列データ14からそれぞれ、特徴行列F1、F2、及び F3 を算出する。そして情報処理装置2000は、これら3つの特徴行列の組を、対象ガスの特徴量として出力する。なお、前述したように、複数の特徴行列の組を出力する代わりに、これら複数の特徴行列を連結した1つの特徴行列Fc を出力してもよい。
 ここで、特性の異なる複数のセンサ10は、1つの筐体に収められてもよいし、それぞれ異なる筐体に収められてもよい。前者の場合、例えば、1つのセンサ筐体の中に材質の異なる複数の官能膜を収納し、各官能膜について検出値が得られるように、センサ10が構成される。
 さらに、図11で説明した方法と、図12で説明した方法を組み合わせてもよい。すなわち、情報処理装置2000は、複数のセンサ10それぞれから、立ち上がりの時系列データ14と立ち下がりの時系列データ14を得て、得られた各時系列データ14について特徴行列 F を算出し、算出した複数の特徴行列の組や、これらを連結した1つの特徴行列を、対象ガスの特徴量としてもよい。
<バイアスを考慮した特徴量の算出>
 センサ10の検出値には、時間に応じた変化を表さないバイアス項が含まれていることがある。この場合、時系列データ14は以下のように表される。なお、ここでは特徴定数として、速度定数βを用いている。
Figure JPOXMLDOC01-appb-M000024
 バイアスは、例えば、センサ10のオフセットがずれていることによって生じる。その他にも例えば、バイアスは、対象ガスとパージガスに共通して含まれている成分の寄与(例えば、大気中の窒素や酸素の寄与)により生じる。
 情報処理装置2000は、時系列データ14からオフセットを除去する機能を有していてもよい。こうすることで、対象ガスの特徴量をより正確に算出することができる。以下、オフセットを考慮して特徴量を算出する方法について説明する。
 算出部2040は、上記式(27)で表される時系列データ14の予測モデルを生成することにより、バイアスを考慮して寄与ベクトルΞを算出する。すなわち、算出部2040は、式(27)で表現された予測モデルについて、パラメータΞ及び b の推定を行う。具体的には、算出部2040は、目的関数(10)、(12)、又は(16)を、Ξのみならず、b についても最適化することにより、Ξと b を推定する。なお、特徴定数として時定数を用いる場合には、式(27)においてβk を1/τk に置き換える。
 例えば、目的関数として式(14)を用いるとする。この場合、算出部2040は、Ξ及び b を以下の最適化問題によって計算する。目的関数に(8)や(10)を用いる場合も同様である。
Figure JPOXMLDOC01-appb-M000025
 上記最適化問題の解Ξ及び b は、以下の式によって算出することができる。
Figure JPOXMLDOC01-appb-M000026
 このように、バイアス b と寄与ベクトルΞの両方を推定することにより、寄与ベクトルからバイアスの効果が取り除かれ、センサ10の検定値にバイアスが含まれている場合にも正確に寄与ベクトルを計算することができる。
 なお、出力部2060は、特徴行列 F に加え、バイアス b や b0 を出力しても良い。バイアスが上記センサのオフセットのずれによるものである場合、b0の値を用いてセンサのオフセットをキャリブレーションすることができる。
 以上、図面を参照して本発明の実施形態について述べたが、これらは本発明の例示であり、上記以外の様々な構成を採用することもできる。
 上記の実施形態の一部又は全部は、以下の付記のようにも記載されうるが、以下には限られない。
1. 対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得部と、
 前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出部と、
 複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力部と、を有し、
 前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、情報処理装置。
2. 前記算出部は、
  前記時系列データの測定期間から複数の部分期間を抽出し、
  前記部分期間ごとに、その部分期間における前記検出値の時間変化率の対数に基づいて、前記特徴定数を算出し、
 前記部分期間は、その中に含まれる検出値の時間変化率の対数が略同一の期間である、1.に記載の情報処理装置。
3. 前記算出部は、
  前記時系列データを用いて、各時刻における検出値及びその時刻における検出値の時間変化率を要素とする時系列のベクトルデータを算出し、
  各前記算出した時系列のベクトルデータそれぞれについて速度ベクトルを算出し、
  前記速度ベクトルの向きに基づいて、前記時系列データの測定期間から、複数の部分期間を抽出し、
  前記部分期間ごとに、その部分期間における前記速度ベクトルの向きに基づいて、前記特徴定数を算出し、
 前記部分期間は、その中に含まれる前記速度ベクトルの向きが略同一の期間である、1.に記載の情報処理装置。
4. 前記算出部は、複数の特徴定数それぞれの寄与値をパラメータとする前記センサの検出値の予測モデルについて、前記取得した時系列データを用いたパラメータ推定を行うことで、各寄与値を算出する、1.乃至3.いずれか一つに記載の情報処理装置。
5. 前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとについて、最小二乗法による最尤推定を行うことで、各前記寄与値を算出する、4.に記載の情報処理装置。
6. 前記最小二乗法における最尤推定において、目的関数に正則化項が含まれている、5.に記載の情報処理装置。
7. 前記算出部は、各前記寄与値の事前分布と前記取得した時系列データとを用いた MAP(Maximum a Posteriori)推定又はベイズ推定により、各前記寄与値を算出する、4.に記載の情報処理装置。
8. 前記事前分布は、多変量正規分布又はガウス過程である、7.に記載の情報処理装置。
9. 前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとの二乗誤差を表す目的関数について、その目的関数の最小値を複数の特徴定数に関して最小化することにより、複数の特徴定数及び複数の寄与値を算出する、4.に記載の情報処理装置。
10. 前記予測モデルには、バイアスを表すパラメータが含まれており、
 前記算出部は、前記予測モデルについて寄与値とバイアスそれぞれを表すパラメータを推定する、4.乃至9.いずれか一つに記載の情報処理装置。
11. 前記時系列データ取得部は、複数の時系列データを取得し、
 前記算出部は、複数の時系列データそれぞれについて、特徴定数の集合と寄与値の集合の組みを算出し、
 前記出力部は、前記算出された特徴定数の集合と寄与値の集合の組みを複数まとめた情報を、前記対象ガスの特徴量として出力する、1.乃至10.いずれか一つに記載の情報処理装置。
12. 前記複数の時系列データは、前記センサを前記対象ガスに曝す際に得られる時系列データと、前記センサから前記対象ガスを取り除く際に得られる時系列データとの双方を含む、11.に記載の情報処理装置。
13. 前記複数の時系列データは、特性の異なる複数の前記センサそれぞれから得られる時系列データを含む、11.に記載の情報処理装置。
14. コンピュータによって実行させる制御方法であって、
 対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得ステップと、
 前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出ステップと、
 複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力ステップと、を有し、
 前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、制御方法。
15. 前記算出ステップにおいて、
  前記時系列データの測定期間から複数の部分期間を抽出し、
  前記部分期間ごとに、その部分期間における前記検出値の時間変化率の対数に基づいて、前記特徴定数を算出し、
 前記部分期間は、その中に含まれる検出値の時間変化率の対数が略同一の期間である、14.に記載の制御方法。
16. 前記算出ステップにおいて、
  前記時系列データを用いて、各時刻における検出値及びその時刻における検出値の時間変化率を要素とする時系列のベクトルデータを算出し、
  各前記算出した時系列のベクトルデータそれぞれについて速度ベクトルを算出し、
  前記速度ベクトルの向きに基づいて、前記時系列データの測定期間から、複数の部分期間を抽出し、
  前記部分期間ごとに、その部分期間における前記速度ベクトルの向きに基づいて、前記特徴定数を算出し、
 前記部分期間は、その中に含まれる前記速度ベクトルの向きが略同一の期間である、14.に記載の制御方法。
17. 前記算出ステップにおいて、複数の特徴定数それぞれの寄与値をパラメータとする前記センサの検出値の予測モデルについて、前記取得した時系列データを用いたパラメータ推定を行うことで、各寄与値を算出する、14.乃至16.いずれか一つに記載の制御方法。
18. 前記算出ステップにおいて、前記予測モデルから得られる時系列データと前記取得した時系列データとについて、最小二乗法による最尤推定を行うことで、各前記寄与値を算出する、17.に記載の制御方法。
19. 前記最小二乗法における最尤推定において、目的関数に正則化項が含まれている、18.に記載の制御方法。
20. 前記算出ステップにおいて、各前記寄与値の事前分布と前記取得した時系列データとを用いた MAP(Maximum a Posteriori)推定又はベイズ推定により、各前記寄与値を算出する、17.に記載の制御方法。
21. 前記事前分布は、多変量正規分布又はガウス過程である、20.に記載の制御方法。
22. 前記算出ステップにおいて、前記予測モデルから得られる時系列データと前記取得した時系列データとの二乗誤差を表す目的関数について、その目的関数の最小値を複数の特徴定数に関して最小化することにより、複数の特徴定数及び複数の寄与値を算出する、17.に記載の制御方法。
23. 前記予測モデルには、バイアスを表すパラメータが含まれており、
 前記算出ステップにおいて、前記予測モデルについて寄与値とバイアスそれぞれを表すパラメータを推定する、17.乃至22.いずれか一つに記載の制御方法。
24. 前記時系列データ取得ステップにおいて、複数の時系列データを取得し、
 前記算出ステップにおいて、複数の時系列データそれぞれについて、特徴定数の集合と寄与値の集合の組みを算出し、
 前記出力ステップにおいて、前記算出された特徴定数の集合と寄与値の集合の組みを複数まとめた情報を、前記対象ガスの特徴量として出力する、14.乃至23.いずれか一つに記載の制御方法。
25. 前記複数の時系列データは、前記センサを前記対象ガスに曝す際に得られる時系列データと、前記センサから前記対象ガスを取り除く際に得られる時系列データとの双方を含む、24.に記載の制御方法。
26. 前記複数の時系列データは、特性の異なる複数の前記センサそれぞれから得られる時系列データを含む、24.に記載の制御方法。
27. 14.乃至26.いずれか一つに記載の制御方法の各ステップをコンピュータに実行させるプログラム。

Claims (27)

  1.  対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得部と、
     前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出部と、
     複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力部と、を有し、
     前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、情報処理装置。
  2.  前記算出部は、
      前記時系列データの測定期間から複数の部分期間を抽出し、
      前記部分期間ごとに、その部分期間における前記検出値の時間変化率の対数に基づいて、前記特徴定数を算出し、
     前記部分期間は、その中に含まれる検出値の時間変化率の対数が略同一の期間である、請求項1に記載の情報処理装置。
  3.  前記算出部は、
      前記時系列データを用いて、各時刻における検出値及びその時刻における検出値の時間変化率を要素とする時系列のベクトルデータを算出し、
      各前記算出した時系列のベクトルデータそれぞれについて速度ベクトルを算出し、
      前記速度ベクトルの向きに基づいて、前記時系列データの測定期間から、複数の部分期間を抽出し、
      前記部分期間ごとに、その部分期間における前記速度ベクトルの向きに基づいて、前記特徴定数を算出し、
     前記部分期間は、その中に含まれる前記速度ベクトルの向きが略同一の期間である、請求項1に記載の情報処理装置。
  4.  前記算出部は、複数の特徴定数それぞれの寄与値をパラメータとする前記センサの検出値の予測モデルについて、前記取得した時系列データを用いたパラメータ推定を行うことで、各寄与値を算出する、請求項1乃至3いずれか一項に記載の情報処理装置。
  5.  前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとについて、最小二乗法による最尤推定を行うことで、各前記寄与値を算出する、請求項4に記載の情報処理装置。
  6.  前記最小二乗法における最尤推定において、目的関数に正則化項が含まれている、請求項5に記載の情報処理装置。
  7.  前記算出部は、各前記寄与値の事前分布と前記取得した時系列データとを用いた MAP(Maximum a Posteriori)推定又はベイズ推定により、各前記寄与値を算出する、請求項4に記載の情報処理装置。
  8.  前記事前分布は、多変量正規分布又はガウス過程である、請求項7に記載の情報処理装置。
  9.  前記算出部は、前記予測モデルから得られる時系列データと前記取得した時系列データとの二乗誤差を表す目的関数について、その目的関数の最小値を複数の特徴定数に関して最小化することにより、複数の特徴定数及び複数の寄与値を算出する、請求項4に記載の情報処理装置。
  10.  前記予測モデルには、バイアスを表すパラメータが含まれており、
     前記算出部は、前記予測モデルについて寄与値とバイアスそれぞれを表すパラメータを推定する、請求項4乃至9いずれか一項に記載の情報処理装置。
  11.  前記時系列データ取得部は、複数の時系列データを取得し、
     前記算出部は、複数の時系列データそれぞれについて、特徴定数の集合と寄与値の集合の組みを算出し、
     前記出力部は、前記算出された特徴定数の集合と寄与値の集合の組みを複数まとめた情報を、前記対象ガスの特徴量として出力する、請求項1乃至10いずれか一項に記載の情報処理装置。
  12.  前記複数の時系列データは、前記センサを前記対象ガスに曝す際に得られる時系列データと、前記センサから前記対象ガスを取り除く際に得られる時系列データとの双方を含む、請求項11に記載の情報処理装置。
  13.  前記複数の時系列データは、特性の異なる複数の前記センサそれぞれから得られる時系列データを含む、請求項11に記載の情報処理装置。
  14.  コンピュータによって実行させる制御方法であって、
     対象ガスに含まれる分子の付着と離脱に応じて検出値が変化するセンサから出力された、検出値の時系列データを取得する時系列データ取得ステップと、
     前記時系列データに対して寄与する複数の特徴定数と、各前記特徴定数の前記時系列データに対する寄与の大きさを表す寄与値とを算出する算出ステップと、
     複数の特徴定数と各特徴定数について算出された寄与値の組み合わせを、前記センサによってセンシングされたガスの特徴量として出力する出力ステップと、を有し、
     前記特徴定数は、前記センサに付着している分子の量の時間変化の大きさに関する時定数又は速度定数である、制御方法。
  15.  前記算出ステップにおいて、
      前記時系列データの測定期間から複数の部分期間を抽出し、
      前記部分期間ごとに、その部分期間における前記検出値の時間変化率の対数に基づいて、前記特徴定数を算出し、
     前記部分期間は、その中に含まれる検出値の時間変化率の対数が略同一の期間である、請求項14に記載の制御方法。
  16.  前記算出ステップにおいて、
      前記時系列データを用いて、各時刻における検出値及びその時刻における検出値の時間変化率を要素とする時系列のベクトルデータを算出し、
      各前記算出した時系列のベクトルデータそれぞれについて速度ベクトルを算出し、
      前記速度ベクトルの向きに基づいて、前記時系列データの測定期間から、複数の部分期間を抽出し、
      前記部分期間ごとに、その部分期間における前記速度ベクトルの向きに基づいて、前記特徴定数を算出し、
     前記部分期間は、その中に含まれる前記速度ベクトルの向きが略同一の期間である、請求項14に記載の制御方法。
  17.  前記算出ステップにおいて、複数の特徴定数それぞれの寄与値をパラメータとする前記センサの検出値の予測モデルについて、前記取得した時系列データを用いたパラメータ推定を行うことで、各寄与値を算出する、請求項14乃至16いずれか一項に記載の制御方法。
  18.  前記算出ステップにおいて、前記予測モデルから得られる時系列データと前記取得した時系列データとについて、最小二乗法による最尤推定を行うことで、各前記寄与値を算出する、請求項17に記載の制御方法。
  19.  前記最小二乗法における最尤推定において、目的関数に正則化項が含まれている、請求項18に記載の制御方法。
  20.  前記算出ステップにおいて、各前記寄与値の事前分布と前記取得した時系列データとを用いた MAP(Maximum a Posteriori)推定又はベイズ推定により、各前記寄与値を算出する、請求項17に記載の制御方法。
  21.  前記事前分布は、多変量正規分布又はガウス過程である、請求項20に記載の制御方法。
  22.  前記算出ステップにおいて、前記予測モデルから得られる時系列データと前記取得した時系列データとの二乗誤差を表す目的関数について、その目的関数の最小値を複数の特徴定数に関して最小化することにより、複数の特徴定数及び複数の寄与値を算出する、請求項17に記載の制御方法。
  23.  前記予測モデルには、バイアスを表すパラメータが含まれており、
     前記算出ステップにおいて、前記予測モデルについて寄与値とバイアスそれぞれを表すパラメータを推定する、請求項17乃至22いずれか一項に記載の制御方法。
  24.  前記時系列データ取得ステップにおいて、複数の時系列データを取得し、
     前記算出ステップにおいて、複数の時系列データそれぞれについて、特徴定数の集合と寄与値の集合の組みを算出し、
     前記出力ステップにおいて、前記算出された特徴定数の集合と寄与値の集合の組みを複数まとめた情報を、前記対象ガスの特徴量として出力する、請求項14乃至23いずれか一項に記載の制御方法。
  25.  前記複数の時系列データは、前記センサを前記対象ガスに曝す際に得られる時系列データと、前記センサから前記対象ガスを取り除く際に得られる時系列データとの双方を含む、請求項24に記載の制御方法。
  26.  前記複数の時系列データは、特性の異なる複数の前記センサそれぞれから得られる時系列データを含む、請求項24に記載の制御方法。
  27.  請求項14乃至26いずれか一項に記載の制御方法の各ステップをコンピュータに実行させるプログラム。
PCT/JP2018/028565 2018-07-31 2018-07-31 情報処理装置、制御方法、及びプログラム WO2020026327A1 (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
JP2020533925A JP7140191B2 (ja) 2018-07-31 2018-07-31 情報処理装置、制御方法、及びプログラム
US17/262,955 US20210311009A1 (en) 2018-07-31 2018-07-31 Information processing apparatus, control method, and non-transitory storage medium
PCT/JP2018/028565 WO2020026327A1 (ja) 2018-07-31 2018-07-31 情報処理装置、制御方法、及びプログラム

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
PCT/JP2018/028565 WO2020026327A1 (ja) 2018-07-31 2018-07-31 情報処理装置、制御方法、及びプログラム

Publications (1)

Publication Number Publication Date
WO2020026327A1 true WO2020026327A1 (ja) 2020-02-06

Family

ID=69232427

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2018/028565 WO2020026327A1 (ja) 2018-07-31 2018-07-31 情報処理装置、制御方法、及びプログラム

Country Status (3)

Country Link
US (1) US20210311009A1 (ja)
JP (1) JP7140191B2 (ja)
WO (1) WO2020026327A1 (ja)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20210311009A1 (en) * 2018-07-31 2021-10-07 Nec Corporation Information processing apparatus, control method, and non-transitory storage medium
WO2023037999A1 (ja) * 2021-09-07 2023-03-16 パナソニックIpマネジメント株式会社 ガス分析方法及びガス分析システム

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH03276046A (ja) * 1990-03-27 1991-12-06 Res Dev Corp Of Japan 多入力型過渡波形解析装置
JPH11264809A (ja) * 1998-03-17 1999-09-28 Shimadzu Corp ガス識別装置
JP2006275606A (ja) * 2005-03-28 2006-10-12 Kyoto Univ ガス検出方法及びガス検出装置
US20090169089A1 (en) * 2005-05-16 2009-07-02 Isis Innovation Limited Cell Analysis
WO2014103560A1 (ja) * 2012-12-25 2014-07-03 インターナショナル・ビジネス・マシーンズ・コーポレーション 分析装置、分析プログラム、分析方法、推定装置、推定プログラム、及び、推定方法。
US20160084808A1 (en) * 2014-09-20 2016-03-24 Commissariat A L'energie Atomique Et Aux Energies Alternatives Method and device for determining a composition of a gas sample processed by means of gas chromatography
WO2018101128A1 (ja) * 2016-11-29 2018-06-07 国立研究開発法人物質・材料研究機構 試料に対応付けられた推定対象値を推定する方法及び装置

Family Cites Families (26)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPWO2004090517A1 (ja) * 2003-04-04 2006-07-06 独立行政法人産業技術総合研究所 蛍光寿命を利用した物質の定量用試薬、方法及び装置
JP5598200B2 (ja) * 2010-09-16 2014-10-01 ソニー株式会社 データ処理装置、データ処理方法、およびプログラム
US9122705B1 (en) * 2011-03-15 2015-09-01 Google Inc. Scoring hash functions
US8732528B1 (en) * 2012-01-06 2014-05-20 Amazon Technologies, Inc. Measuring test effects using adjusted outlier data
US10309924B2 (en) * 2013-06-07 2019-06-04 Cornell University Floating gate based sensor apparatus and related floating gate based sensor applications
EP3062103A1 (en) * 2015-02-27 2016-08-31 Alpha M.O.S. Portable fluid sensory device with learning capabilities
US10127694B2 (en) * 2016-11-18 2018-11-13 Adobe Systems Incorporated Enhanced triplet embedding and triplet creation for high-dimensional data visualizations
JP7069581B2 (ja) * 2017-07-07 2022-05-18 富士通株式会社 分析方法、分析装置及び分析プログラム
JP6950504B2 (ja) * 2017-12-08 2021-10-13 富士通株式会社 異常候補抽出プログラム、異常候補抽出方法および異常候補抽出装置
US20210232957A1 (en) * 2018-06-07 2021-07-29 Nec Corporation Relationship analysis device, relationship analysis method, and recording medium
JP7017711B2 (ja) * 2018-06-07 2022-02-09 日本電気株式会社 関係性分析装置、関係性分析方法およびプログラム
JP7074194B2 (ja) * 2018-07-31 2022-05-24 日本電気株式会社 情報処理装置、制御方法、及びプログラム
WO2020026327A1 (ja) * 2018-07-31 2020-02-06 日本電気株式会社 情報処理装置、制御方法、及びプログラム
JP7101084B2 (ja) * 2018-08-29 2022-07-14 株式会社東芝 情報処理装置、情報処理システム及び情報処理方法
JP7063389B2 (ja) * 2018-09-27 2022-05-09 日本電気株式会社 処理装置、処理方法、およびプログラム
WO2020100285A1 (ja) * 2018-11-16 2020-05-22 日本電気株式会社 情報処理装置、制御方法、及びプログラム
CN113424199A (zh) * 2019-01-23 2021-09-21 谷歌有限责任公司 用于神经网络的复合模型缩放
WO2020202338A1 (ja) * 2019-03-29 2020-10-08 日本電気株式会社 情報処理装置、情報処理方法、およびプログラム
US20220309397A1 (en) * 2019-06-19 2022-09-29 Nec Corporation Prediction model re-learning device, prediction model re-learning method, and program recording medium
US11468364B2 (en) * 2019-09-09 2022-10-11 Humana Inc. Determining impact of features on individual prediction of machine learning based models
JP7384059B2 (ja) * 2020-02-06 2023-11-21 富士通株式会社 検知プログラム、検知方法及び検知装置
US20230214668A1 (en) * 2020-09-10 2023-07-06 Aizoth Inc. Hyperparameter adjustment device, non-transitory recording medium in which hyperparameter adjustment program is recorded, and hyperparameter adjustment program
CN112418395B (zh) * 2020-11-17 2022-08-26 吉林大学 一种基于生成对抗网络的气体传感器阵列漂移补偿方法
CN112464999B (zh) * 2020-11-19 2021-10-19 深圳大学 一种混合气体种类与浓度识别方法、装置及存储介质
US20220172086A1 (en) * 2020-11-30 2022-06-02 Nice Ltd. System and method for providing unsupervised model health monitoring
CN117897611A (zh) * 2021-09-07 2024-04-16 松下知识产权经营株式会社 气体分析方法以及气体分析系统

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH03276046A (ja) * 1990-03-27 1991-12-06 Res Dev Corp Of Japan 多入力型過渡波形解析装置
JPH11264809A (ja) * 1998-03-17 1999-09-28 Shimadzu Corp ガス識別装置
JP2006275606A (ja) * 2005-03-28 2006-10-12 Kyoto Univ ガス検出方法及びガス検出装置
US20090169089A1 (en) * 2005-05-16 2009-07-02 Isis Innovation Limited Cell Analysis
WO2014103560A1 (ja) * 2012-12-25 2014-07-03 インターナショナル・ビジネス・マシーンズ・コーポレーション 分析装置、分析プログラム、分析方法、推定装置、推定プログラム、及び、推定方法。
US20160084808A1 (en) * 2014-09-20 2016-03-24 Commissariat A L'energie Atomique Et Aux Energies Alternatives Method and device for determining a composition of a gas sample processed by means of gas chromatography
WO2018101128A1 (ja) * 2016-11-29 2018-06-07 国立研究開発法人物質・材料研究機構 試料に対応付けられた推定対象値を推定する方法及び装置

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
SETO, SHUICHI ET AL.: "Chemical Analysis for Transient Response of Semiconductor Sensor by Autoregressive Model", IEEJ TRANSACTIONS ON SENSORS AND MICROMACHINES, vol. 125, no. 3, 2005, pages 129 - 134 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20210311009A1 (en) * 2018-07-31 2021-10-07 Nec Corporation Information processing apparatus, control method, and non-transitory storage medium
WO2023037999A1 (ja) * 2021-09-07 2023-03-16 パナソニックIpマネジメント株式会社 ガス分析方法及びガス分析システム

Also Published As

Publication number Publication date
US20210311009A1 (en) 2021-10-07
JP7140191B2 (ja) 2022-09-21
JPWO2020026327A1 (ja) 2021-08-02

Similar Documents

Publication Publication Date Title
Quinn et al. Understanding sequencing data as compositions: an outlook and review
JP7063389B2 (ja) 処理装置、処理方法、およびプログラム
TW201224431A (en) Monitoring, detecting and quantifying chemical compounds in a sample
EP2035803A1 (en) Sensor system for estimating varying field
US20210232957A1 (en) Relationship analysis device, relationship analysis method, and recording medium
WO2020026327A1 (ja) 情報処理装置、制御方法、及びプログラム
JP7276450B2 (ja) 予測モデル再学習装置、予測モデル再学習方法及びプログラム
JP7074194B2 (ja) 情報処理装置、制御方法、及びプログラム
JP7127697B2 (ja) 情報処理装置、制御方法、及びプログラム
JP2019016193A (ja) 分析方法、分析装置及び分析プログラム
WO2019026193A1 (ja) 情報処理装置、情報処理システム、情報処理方法、及び、記録媒体
JP7099623B2 (ja) 情報処理装置、情報処理方法、およびプログラム
JP7056747B2 (ja) 情報処理装置、処理装置、情報処理方法、処理方法、決定方法、およびプログラム
JP2021128423A (ja) 異常検知装置及び異常検知プログラム
WO2020026328A1 (ja) 情報処理装置、制御方法、及びプログラム
JP7006799B2 (ja) 情報処理装置、センサ動作最適化方法、およびプログラム
EP3276343B1 (en) Method and device for characterising an analyte
WO2020065890A1 (ja) 情報処理装置、情報処理方法、およびプログラム
JP7143894B2 (ja) 情報処理装置、伝達関数生成方法、およびプログラム
JP2018151913A (ja) 情報処理システム、情報処理方法、及びプログラム
JP2005267474A (ja) マハラノビス距離を利用した異常原因診断方法及びプログラム
Waller Using surrogate data for nonlinear identification: A case study
Stark Analysis of time series
JP2012093423A (ja) パターン認識方法

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

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2020533925

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

Country of ref document: EP

Kind code of ref document: A1