WO2021024351A1 - 異常推定装置、異常検出装置、異常検出方法、異常検出プログラム - Google Patents

異常推定装置、異常検出装置、異常検出方法、異常検出プログラム Download PDF

Info

Publication number
WO2021024351A1
WO2021024351A1 PCT/JP2019/030717 JP2019030717W WO2021024351A1 WO 2021024351 A1 WO2021024351 A1 WO 2021024351A1 JP 2019030717 W JP2019030717 W JP 2019030717W WO 2021024351 A1 WO2021024351 A1 WO 2021024351A1
Authority
WO
WIPO (PCT)
Prior art keywords
section
data
abnormality
series data
time series
Prior art date
Application number
PCT/JP2019/030717
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 JP2021538566A priority Critical patent/JP7248127B2/ja
Priority to PCT/JP2019/030717 priority patent/WO2021024351A1/ja
Priority to US17/632,527 priority patent/US12090940B2/en
Publication of WO2021024351A1 publication Critical patent/WO2021024351A1/ja

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01P—MEASURING LINEAR OR ANGULAR SPEED, ACCELERATION, DECELERATION, OR SHOCK; INDICATING PRESENCE, ABSENCE, OR DIRECTION, OF MOVEMENT
    • G01P15/00—Measuring acceleration; Measuring deceleration; Measuring shock, i.e. sudden change of acceleration
    • G01P15/02—Measuring acceleration; Measuring deceleration; Measuring shock, i.e. sudden change of acceleration by making use of inertia forces using solid seismic masses
    • G01P15/08—Measuring acceleration; Measuring deceleration; Measuring shock, i.e. sudden change of acceleration by making use of inertia forces using solid seismic masses with conversion into electric or magnetic values
    • G01P15/0891—Measuring acceleration; Measuring deceleration; Measuring shock, i.e. sudden change of acceleration by making use of inertia forces using solid seismic masses with conversion into electric or magnetic values with indication of predetermined acceleration values
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B60—VEHICLES IN GENERAL
    • B60R—VEHICLES, VEHICLE FITTINGS, OR VEHICLE PARTS, NOT OTHERWISE PROVIDED FOR
    • B60R21/00—Arrangements or fittings on vehicles for protecting or preventing injuries to occupants or pedestrians in case of accidents or other traffic risks
    • B60R21/01—Electrical circuits for triggering passive safety arrangements, e.g. airbags, safety belt tighteners, in case of vehicle accidents or impending vehicle accidents
    • B60R21/013—Electrical circuits for triggering passive safety arrangements, e.g. airbags, safety belt tighteners, in case of vehicle accidents or impending vehicle accidents including means for detecting collisions, impending collisions or roll-over
    • B60R21/0132—Electrical circuits for triggering passive safety arrangements, e.g. airbags, safety belt tighteners, in case of vehicle accidents or impending vehicle accidents including means for detecting collisions, impending collisions or roll-over responsive to vehicle motion parameters, e.g. to vehicle longitudinal or transversal deceleration or speed value
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B60—VEHICLES IN GENERAL
    • B60W—CONJOINT CONTROL OF VEHICLE SUB-UNITS OF DIFFERENT TYPE OR DIFFERENT FUNCTION; CONTROL SYSTEMS SPECIALLY ADAPTED FOR HYBRID VEHICLES; ROAD VEHICLE DRIVE CONTROL SYSTEMS FOR PURPOSES NOT RELATED TO THE CONTROL OF A PARTICULAR SUB-UNIT
    • B60W30/00—Purposes of road vehicle drive control systems not related to the control of a particular sub-unit, e.g. of systems using conjoint control of vehicle sub-units
    • B60W30/08—Active safety systems predicting or avoiding probable or impending collision or attempting to minimise its consequences
    • G—PHYSICS
    • G07—CHECKING-DEVICES
    • G07C—TIME OR ATTENDANCE REGISTERS; REGISTERING OR INDICATING THE WORKING OF MACHINES; GENERATING RANDOM NUMBERS; VOTING OR LOTTERY APPARATUS; ARRANGEMENTS, SYSTEMS OR APPARATUS FOR CHECKING NOT PROVIDED FOR ELSEWHERE
    • G07C5/00—Registering or indicating the working of vehicles
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B60—VEHICLES IN GENERAL
    • B60R—VEHICLES, VEHICLE FITTINGS, OR VEHICLE PARTS, NOT OTHERWISE PROVIDED FOR
    • B60R21/00—Arrangements or fittings on vehicles for protecting or preventing injuries to occupants or pedestrians in case of accidents or other traffic risks
    • B60R21/01—Electrical circuits for triggering passive safety arrangements, e.g. airbags, safety belt tighteners, in case of vehicle accidents or impending vehicle accidents
    • B60R21/013—Electrical circuits for triggering passive safety arrangements, e.g. airbags, safety belt tighteners, in case of vehicle accidents or impending vehicle accidents including means for detecting collisions, impending collisions or roll-over
    • B60R2021/01302—Electrical circuits for triggering passive safety arrangements, e.g. airbags, safety belt tighteners, in case of vehicle accidents or impending vehicle accidents including means for detecting collisions, impending collisions or roll-over monitoring vehicle body vibrations or noise
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B60—VEHICLES IN GENERAL
    • B60W—CONJOINT CONTROL OF VEHICLE SUB-UNITS OF DIFFERENT TYPE OR DIFFERENT FUNCTION; CONTROL SYSTEMS SPECIALLY ADAPTED FOR HYBRID VEHICLES; ROAD VEHICLE DRIVE CONTROL SYSTEMS FOR PURPOSES NOT RELATED TO THE CONTROL OF A PARTICULAR SUB-UNIT
    • B60W2400/00—Indexing codes relating to detected, measured or calculated conditions or factors
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01H—MEASUREMENT OF MECHANICAL VIBRATIONS OR ULTRASONIC, SONIC OR INFRASONIC WAVES
    • G01H1/00—Measuring characteristics of vibrations in solids by using direct conduction to the detector
    • G01H1/003—Measuring characteristics of vibrations in solids by using direct conduction to the detector of rotating machines
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01H—MEASUREMENT OF MECHANICAL VIBRATIONS OR ULTRASONIC, SONIC OR INFRASONIC WAVES
    • G01H17/00—Measuring mechanical vibrations or ultrasonic, sonic or infrasonic waves, not provided for in the preceding groups

Definitions

  • the present invention relates to an anomaly estimation device, anomaly detection device, anomaly detection method, and anomaly detection program in time-series data in which a feature pattern can change with the passage of time.
  • An example of time-series data in which the feature pattern can change with the passage of time is data related to vibration acquired by a sensor mounted on an automobile.
  • a typical example of an in-vehicle sensor is a drive recorder, and when it detects an acceleration that exceeds a preset threshold value (for example, 0.2G), it determines that there was a collision and detects an abnormality. ing.
  • a preset threshold value for example, 0.2G
  • the present invention handles time-series data related to vibration acquired by an in-vehicle sensor as typical time-series data.
  • the time-series data whose feature pattern can change with the passage of time is not limited to the time-series data related to vibration acquired by the in-vehicle sensor.
  • Techniques for identifying an abnormality from general time-series data are called abnormality detection and change detection, and several methods have been proposed in Non-Patent Document 1 and the like.
  • Non-Patent Document 2 discloses a technique related to dynamic mode analysis.
  • Non-Patent Documents 1 and 2 has been changed to a technology for detecting anomalies from time-series data such as time-series data related to automobile vibration, in which the feature pattern may change with the passage of time. It is unknown if it is available.
  • a first object of the present invention is to provide an abnormality estimation device capable of limiting the images that must be confirmed when investigating the cause of a scratch on a vehicle driven by a plurality of users. ..
  • a second object of the present invention is to provide an abnormality detection device capable of acquiring an abnormality degree required for abnormality detection even if the feature pattern of the time series data changes with the passage of time.
  • the abnormality estimation device of the present invention is an abnormality estimation device for estimating the time when an abnormality occurs in a vehicle driven by a plurality of users.
  • the abnormality estimation device of the present invention includes a sensor, a prediction unit, an abnormality degree acquisition unit, and an abnormality determination unit.
  • the sensor is installed in the car and acquires time series data related to vibration.
  • the prediction unit uses the data of the first section, which is the data of the first section, which is a part of the time series data, to the data of the second section, which is a part of the time series data after the first section. Predict and obtain the predicted second section data.
  • the abnormality degree acquisition unit acquires the abnormality degree in the second section based on the difference between the predicted second section data and the second section data which is the actual data of the second section in the time series data.
  • the abnormality determination unit determines whether or not an abnormality has occurred in the second section based on the acquired abnormality degree.
  • the abnormality detection device of the present invention is an abnormality detection device that detects the degree of abnormality in time series data.
  • the abnormality detection device of the present invention includes a first acquisition unit, a prediction unit, and an abnormality degree acquisition unit.
  • the first acquisition unit acquires the dynamic feature pattern of the first section from the first section, which is a part of the time series data.
  • the prediction unit predicts the data of the second section, which is a part of the time series data after the first section, using the feature pattern, and obtains the prediction second section data.
  • the abnormality degree acquisition unit acquires the abnormality degree based on the difference between the predicted second section data and the actual data of the second section in the time series data.
  • the abnormality estimation device of the present invention at least the range where there is no abnormality can be excluded, so that it is possible to limit the images to be confirmed when investigating the cause of the scratches on the automobile driven by a plurality of users.
  • the abnormality detection device of the present invention the predicted second section data is obtained based on the dynamic feature pattern, so that even if the feature pattern of the time series data changes with the passage of time, it is necessary for abnormality detection. You can get the degree of abnormality.
  • FIG. 1 is an example of the probability distribution at the time of normal and the time of abnormality.
  • FIG. 2 shows the method of candidate 2.
  • This is a method using the singular spectrum transformation shown in Non-Patent Document 1.
  • the degree of difference (probability density ratio) between the past and present probability distributions is evaluated using the subspace method.
  • characteristic patterns are extracted by singular value decomposition of time series data of past and present sections, and the subspace where the past feature pattern is stretched and the subspace where the current feature pattern is stretched.
  • the degree of change is calculated from the size of the overlap of.
  • the current section means a section for which the degree of change is to be obtained. Specifically, the degree of change is calculated as follows.
  • Figure 3 shows the procedure for creating a data matrix.
  • K, M, N are integers of 2 or more, n is an integer of 1 or more and N or less, L is an integer of 1 or more, k is an integer of 1 or more and K or less, i, r is an integer of 1 or more, t is 1 or more and M + N An integer less than or equal to -1, [] is a symbol indicating a row vector in which numerical values are arranged in the row direction or a matrix in which column vectors are arranged in the row direction, T is a symbol indicating translocation, and x (i) is time series data.
  • the time series data x (1), ..., X (K (M + N-1)) is divided into M + N-1 sections, and the sections D 1 , D 2 , ..., D k , ... , DK is generated.
  • D k is a section whose section data is x (1+ (k-1) (M + N-1)), ..., X (k (M + N-1)).
  • D k is
  • D k [x (1 + (k-1) (M + N-1)), ..., x (k (M + N-1))] T It may be expressed as a column vector having M + N-1 elements as in.
  • the section data D k is divided into a partial time series having a window length (length of a partial time series) M and a shift width of 1, and N partial time series are created.
  • the nth partial time series x k, n of the kth interval D k is M time series data x (n + (k-1) (M + N-1)), ..., X (n + (k-1)). ) (M + N-1) + M-1) as an element
  • x k, n [x (n + (k-1) (M + N-1)), ..., x (n + (k-1) (M + N-1) + M-1)] T It is represented by a column vector as in. In addition, because there is a limitation in the notation, the column vector is expressed as above, but this notation is
  • N sub-time series x k, 1 , ..., X k, N in the interval D k are put together to form a data matrix X k .
  • x k, 1, x k, n, x k, N indicates that it is a column vector.
  • FIG. 4 shows a method for calculating the degree of change in Candidate 2.
  • L is called a lag and represents the mutual position of the past section (first section) and the current section (second section).
  • Data matrix X k-L of the first section and a second section performs singular value decomposition on X k, take out the singular value r Largest left singular vectors.
  • the matrix of the singular vectors uk -L, 1 , uk -L, 2 , ..., uk -L, r extracted from the data matrix X kl of the first interval is defined below.
  • U k-L, r [uk -L, 1 , u k-L, 2 , ..., u k-L, r ]
  • uk -L, 1 , uk -L, 2 , ..., uk -L, r are column vectors.
  • Q k, r [q k, 1 , q k, 2 , ..., q k, r ]
  • q k, 1 , q k, 2 , ..., Q k, r are column vectors. Note that r may be set to a different value in the first section and the second section.
  • c (k) 1- ⁇ U k-L, r T Q k, r ⁇ 2
  • ⁇ ⁇ 2 represents matrix 2-norm
  • ⁇ U k-L, r T Q k, r ⁇ 2 is U k-L, r T Q k, equal to the largest singular value of r.
  • Fig. 5 shows the change detection result using singular value spectrum conversion for the time series data observed in the actual automobile.
  • the vertical axis of FIG. 5 (A) indicates the amplitude of the input time series data
  • the vertical axis of FIG. 5 (B) indicates the degree of change
  • the horizontal axis of both indicates time.
  • the sampling frequency of the time-series data used in the experiment is 16 kHz
  • the range indicated by the arrow in FIG. 5 (A) is the section in which the abnormality should be detected.
  • the extraction number r of the singular vector was set as the upper number until the cumulative contribution rate of the singular value exceeded 80%.
  • the threshold value of the degree of change is set to 0.4, and the points where the degree of change exceeds the threshold value are marked with a triangle to emphasize the points judged to be abnormal. From FIG. 5, it can be seen that there are some points judged to be abnormal other than the section in which the abnormality should be detected. Further, when FIG. 5 (A) is confirmed, there are large vibrations around 43 seconds, 81 seconds, and 99 seconds. From the result of FIG. 5B, it can be seen that even if the threshold value is changed, it is not possible to detect only these three locations.
  • Non-Patent Document 1 a feature pattern is extracted by singular value decomposition in the section of interest. This corresponds to decomposing the eigenvalues of the variance matrix of the data matrix of the section of interest. That is, since the direction vector in which the variance becomes large is extracted by focusing only on the section of interest, only the static feature pattern assuming stationarity in the section of interest can be extracted. However, in an environment where the situation changes from moment to moment, stationarity cannot be assumed, so it is considered that features cannot be extracted well, the error in the degree of change becomes large, and false detections increase.
  • FIG. 6 shows an example of the functional configuration of the abnormality estimation device
  • FIG. 7 shows an example of the processing flow of the abnormality estimation device.
  • FIG. 8 shows the method of the present invention.
  • the abnormality estimation device 10 is a device for estimating the time when an abnormality occurs in an automobile driven by a plurality of users.
  • the abnormality estimation device 10 includes a sensor 210, a first acquisition unit 110, a prediction unit 120, an abnormality degree acquisition unit 130, and an abnormality determination unit 220.
  • the sensor 210 is installed in an automobile and acquires time series data x (i) related to vibration (S210).
  • the first section past section
  • the second section current section
  • the first acquisition unit 110 acquires a dynamic characteristic pattern of D k-L (S110).
  • the time series data x (i) is divided into M + N-1 sections, and the sections D 1 , D 2 , ..., D k having the time series data as elements. , ..., and D K.
  • the klth section DkL is the first section (past section), and the kth section Dk is the second section (current section).
  • the first acquisition unit 110 when the n-th part of the section D k-L sequence x k-L, n is
  • a subtime series x k, 1 , ..., X k, N which is a column vector, is generated so as to be.
  • the matrix XkL, 1 and the matrix XkL, 2 are expressed as described above in relation to the notation, but these notations are expressed as above.
  • the first acquisition unit 110 using the X k-L, pseudo 1 was low rank approximation calculated on the basis of the singular value decomposition as a result of the inverse matrix X k-L, 1 +,
  • A ⁇ k-L ⁇ k-L ⁇ k-L -1 Request eigenmode [Phi k-L and eigenvalues lambda k-L to be a (S110).
  • X k-L, 1 + V k-L, 1, r ⁇ k-L, 1, r -1 U k-L, 1, r T
  • the linear operator A is subjected to the singular value decomposition of XkL , 1 and the low-rank approximation.
  • A X kL, 2 V k-L, 1, r ⁇ k-L, 1, r -1 U k-L, 1, r T Can be obtained as follows. Calculation cost can be reduced by low-rank approximation in this way.
  • ⁇ k-L X k-L, 2 V k-L, 1, r ⁇ k-L, 1, r -1 W You can ask for it like this.
  • Eigenmode ⁇ k-L and eigenvalues eigenvalues ⁇ k-L, which was determined in this manner is,
  • Prediction unit 120 when based on the first section data which is data of the first section D k-L is a partial section of the sequence data, a part of the time series data after the first section D k-L The data of the second section Dk , which is an interval, is predicted, and the predicted second section data is obtained (S120).
  • the data of the second section predicted using the first partial time series x k, 1 and characteristic pattern of the second section (eigenmode [Phi k-L and eigenvalues eigenvalues lambda k-L), the predicted The second section data ⁇ x (1+ (k-1) (M + N-1)), ..., ⁇ X (N + (k-1) (M + N-1) + M-1) are obtained.
  • eigenmode [Phi k-L and eigenvalues eigenvalues lambda k-L are acquired.
  • ⁇ X k, n ⁇ k-L ⁇ k-L n-1 ⁇ k-L -1 x k, 1 Ask like.
  • the prediction unit 120 first
  • the prediction unit 120 the partial time series x k, 1 and predicted partial time series ⁇ x k, 2, ..., ⁇ x k, the first element of the N-1 respectively from the first predicted second section data Data up to the N-1th ⁇ x (1+ (k-1) (M + N-1)), ..., ⁇ X (N-1 + (k-1) (M + N-1)), and the predicted partial series ⁇ x Predict all elements of k and N Data from Nth to M + N-1th of the second interval data ⁇ x (N + (k-1) (M + N-1)), ..., ⁇ x (N + (k-1)) It may be (M + N-1) + M-1).
  • the abnormality degree acquisition unit 130 acquires the abnormality degree in the second section based on the difference between the predicted second section data and the second section data which is the actual data of the second section in the time series data (S130). That is, the abnormality degree acquisition unit 130 uses the prediction second section data ⁇ x (1+ (k-1) (M + N-1)), ..., ⁇ X (N + (k-1) (M + N-1) + M-1). And, from the actual data x (1+ (k-1) (M + N-1)), ..., x (N + (k-1) (M + N-1) + M-1) of the second section in the time series data, it is abnormal. Degree c (k)
  • the abnormality determination unit 220 determines whether an abnormality has occurred in the second section based on the acquired abnormality degree (S220). For example, a threshold value for the degree of abnormality may be set, and it may be determined that an abnormality has occurred when there is a difference exceeding the threshold value.
  • FIG. 9 shows the result of estimating the abnormality by the procedure of the abnormality estimation device 10 of the present invention with respect to the time series data observed in the actual automobile.
  • the vertical axis of FIG. 9 (A) indicates the amplitude of the input time series data
  • the vertical axis of FIG. 9 (B) indicates the degree of abnormality
  • the horizontal axis of both indicates time.
  • the data in FIG. 9 (A) is the same as in FIG. 5 (A). Similar to the experiment shown in FIG. 5, the sampling frequency of the time series data used in the experiment is 16 kHz, and the range indicated by the arrow in FIG. 9 (A) is the section in which the abnormality should be detected.
  • the rank number r was set as a high-ranking number until the cumulative contribution rate of the singular value exceeded 80%.
  • the threshold value of the degree of abnormality is set to 0.3, and the points where the degree of abnormality exceeds the threshold value are marked with a triangle to emphasize the points judged to be abnormal.
  • the setting of this threshold value corresponds to the above-mentioned abnormality determination unit 220.
  • FIG. 9B it can be seen that there is a point where the abnormality is intensively determined in the section where the abnormality should be detected.
  • FIG. 9B it can be seen that there is a point where the abnormality is intensively determined in the section where the abnormality should be detected.
  • the abnormality estimation device of the present invention it is possible to exclude the range where there is no abnormality, so that it is possible to limit the images to be confirmed when investigating the cause of the scratches on the automobile driven by a plurality of users. ..
  • the first embodiment is limited to automobiles driven by a plurality of users such as a share car and a rental car.
  • the present invention can be widely applied to time series data in which the feature pattern may change with the passage of time. That is, the data acquired by the sensor does not have to be limited to the data related to vibration. Therefore, in the second embodiment, an abnormality detection device that inputs some time-series data and outputs an abnormality degree that is easy to detect an abnormality will be described.
  • FIG. 10 shows an example of the functional configuration of the abnormality detection device of the present invention
  • FIG. 11 shows an example of the processing flow of the abnormality detection device of the present invention
  • the abnormality detection device 20 is a device for detecting the degree of abnormality in time series data.
  • the abnormality detection device 20 includes a first acquisition unit 110, a prediction unit 120, and an abnormality degree acquisition unit 130.
  • the abnormality estimation device 10 the sensor 210, and the abnormality determination unit 220 are shown by dotted lines in FIG.
  • Time-series data x (i) is input to the abnormality detection device 20 (S211).
  • the first acquisition unit 110 when the first section (last section) is a partial section of the series data to obtain a dynamic characteristic pattern of the first section from D k-L (S110).
  • Prediction unit 120 when based on the first section data which is data of the first section D k-L is a partial section of the sequence data, a part of the time series data after the first section D k-L
  • the data of the second section Dk which is an interval, is predicted, and the predicted second section data is obtained (S120).
  • the abnormality degree acquisition unit 130 acquires the abnormality degree in the second section based on the difference between the predicted second section data and the second section data which is the actual data of the second section in the time series data (S130). Steps S110, S120, and S130 are the same as in the first embodiment.
  • the eigenvalue representing the attenuation (increase) rate is not considered, so that they are the same. Even if there is a natural mode with a frequency, if it corresponds to a characteristic pattern that attenuates in the past section and a characteristic pattern that increases in the current section, it is indistinguishable by comparing only the natural mode, so it can be used for abnormality detection. The degree of change that can be used could not be obtained.
  • the dynamic feature pattern of the section of interest can be estimated only from the time series data, the accuracy of estimating the degree of change of the time series data obtained from the system in which the situation changes dynamically is improved. , The detection accuracy of the changed part is improved. Since the feature pattern extracted in the observation section is a signal corresponding to the natural mode (natural frequency), the data can be easily interpreted. Therefore, according to the abnormality detection device of the present invention, the prediction second interval data is obtained based on the dynamic feature pattern, so that even if the feature pattern of the time series data changes with the passage of time, the abnormality can be detected. You can get the required degree of abnormality.
  • the program that describes this processing content can be recorded on a computer-readable recording medium.
  • the computer-readable recording medium may be, for example, a magnetic recording device, an optical disk, a photomagnetic recording medium, a semiconductor memory, or the like.
  • this program is carried out, for example, by selling, transferring, renting, etc., a portable recording medium such as a DVD or CD-ROM on which the program is recorded.
  • the program may be stored in the storage device of the server computer, and the program may be distributed by transferring the program from the server computer to another computer via a network.
  • a computer that executes such a program first stores, for example, a program recorded on a portable recording medium or a program transferred from a server computer in its own storage device. Then, at the time of executing the process, the computer reads the program stored in its own recording medium and executes the process according to the read program. Further, as another execution form of this program, a computer may read the program directly from a portable recording medium and execute processing according to the program, and further, the program is transferred from the server computer to this computer. It is also possible to execute the process according to the received program one by one each time.
  • ASP Application Service Provider
  • the program in this embodiment includes information used for processing by a computer and equivalent to the program (data that is not a direct command to the computer but has a property of defining the processing of the computer, etc.).
  • the present device is configured by executing a predetermined program on the computer, but at least a part of these processing contents may be realized by hardware.
  • Abnormality estimation device 20
  • Abnormality detection device 110
  • First acquisition unit 120
  • Prediction unit 130
  • Abnormality acquisition unit 210
  • Sensor 220 Abnormality determination unit

Landscapes

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

Abstract

時間の経過に応じて時系列データの特徴パターンが変化したとしても、異常検出に必要な異常度を取得できる異常検出装置を提供する。本発明の異常検出装置は、時系列データにおける異常度を検出する装置である。異常検出装置は、第一取得部、予測部、異常度取得部を備える。第一取得部は、時系列データの一部区間である第一区間から第一区間の動的な特徴パターンを取得する。予測部は、第一区間よりも後の時系列データの一部区間である第二区間のデータを、特徴パターンを用いて予測し、予測第二区間データを得る。異常度取得部は、予測第二区間データと、時系列データにおける第二区間の実際のデータと、の差に基づき異常度を取得する。

Description

異常推定装置、異常検出装置、異常検出方法、異常検出プログラム
　本発明は、時間の経過に応じて特徴パターンが変化しうる時系列データにおける異常推定装置、異常検出装置、異常検出方法、異常検出プログラムに関する。
　時間の経過に応じて特徴パターンが変化しうる時系列データの例としては、自動車に搭載されたセンサで取得する振動に関するデータなどがある。車載のセンサの代表的な例がドライブレコーダであり、あらかじめ設定されたしきい値（例えば、０．２Ｇなど）を超える加速度を検出すると、衝突があったと判断し、異常を検知する仕組みとなっている。
　本発明では、代表的な時系列データとして、車載のセンサで取得された振動に関する時系列データを扱う。ただし、時間の経過に応じて特徴パターンが変化しうる時系列データであれば、車載のセンサで取得された振動に関する時系列データに限定するものではない。一般的な時系列データから異常があった箇所を特定する技術は、異常検知、変化検知と呼ばれ、非特許文献１などにいくつかの方法が提案されている。また、非特許文献２には、動的モード分析に関する技術が示されている。
井出剛，杉本将，"異常検知と変化検知，"機械学習プロフェッショナルシリーズ，講談社，2015. J. H. Tu, C. W. Rowley, D. M. Luchtenburg, S. L. Brunton, and J. N. Kutz, "On Dynamic Mode Decomposition: Theory and Applications," JCD, vol. 1, no. 2, pp.391-421, Dec. 2014.
　しかしながら、ドライブレコーダの場合、しきい値を超えないような弱い衝突イベント（例えば、車体側面を壁に擦った、車体底面を縁石に擦ったなど）を異常として検知することは困難である。したがって、ドライブレコーダの映像の全てを確認しなければ、車体に傷がついた可能性がある日時を特定できない。したがって、複数の利用者が運転する自動車（例えば、シェアカー、レンタカー、社用車など）の場合、自動車の管理者が車体の傷を見つけても、誰が運転していたかを確認するには、傷が付いたときの映像を見つけるまで、連続的に映像を確認する必要がある。また、非特許文献１，２に示された技術を、自動車の振動に関する時系列データのような時間の経過に応じて特徴パターンが変化する可能性がある時系列データから異常を検知する技術に利用できるのかは不明である。
　本発明は、複数の利用者が運転する自動車に付いた傷の原因を調査する際に、確認しなければならない映像を限定することができる異常推定装置を提供することを第１の目的とする。また、本発明は、時間の経過に応じて時系列データの特徴パターンが変化したとしても、異常検出に必要な異常度を取得できる異常検出装置を提供することを第２の目的とする。
　本発明の異常推定装置は、複数の利用者が運転する自動車において異常が発生した時刻を推定するための異常推定装置である。本発明の異常推定装置は、センサ、予測部、異常度取得部、異常判定部を備える。センサは、自動車に設置され、振動に関する時系列データを取得する。予測部は、時系列データの一部区間である第一区間のデータである第一区間データに基づいて、第一区間よりも後の時系列データの一部区間である第二区間のデータを予測し、予測第二区間データを得る。異常度取得部は、予測第二区間データと、時系列データにおける第二区間の実際のデータである第二区間データと、の差に基づき第二区間における異常度を取得する。異常判定部は、取得された異常度に基づき第二区間において異常が発生しているかを判定する。
　本発明の異常検出装置は、時系列データにおける異常度を検出する異常検出装置である。本発明の異常検出装置は、第一取得部、予測部、異常度取得部を備える。第一取得部は、時系列データの一部区間である第一区間から第一区間の動的な特徴パターンを取得する。予測部は、第一区間よりも後の時系列データの一部区間である第二区間のデータを、特徴パターンを用いて予測し、予測第二区間データを得る。異常度取得部は、予測第二区間データと、時系列データにおける第二区間の実際のデータと、の差に基づき異常度を取得する。
　本発明の異常推定装置によれば、少なくとも異常がない範囲を除外できるので、複数の利用者が運転する自動車に付いた傷の原因を調査する際に、確認しなければならない映像を限定できる。本発明の異常検出装置によれば、動的な特徴パターンに基づいて予測第二区間データを得るので、時間の経過に応じて時系列データの特徴パターンが変化したとしても、異常検出に必要な異常度を取得できる。
正常時と異常時の確率分布の例を示す図。 特異スペクトル変換を用いる手法の概要を示す図。 データ行列の作成手順を示す図。 特異スペクトル変換を用いた変化度の算出手法を示す図。 実際の自動車内で観測された時系列データに対して特異値スペクトル変換を用いた変化検知結果を示す図。 異常推定装置の機能構成例を示す図。 異常推定装置の処理フローの例を示す図。 本発明の手法を示す図。 実際の自動車内で観測された時系列データに対して、本発明の異常推定装置の手順で異常を推定した結果を示す図。 異常検出装置の機能構成例を示す図。 異常検出装置の処理フローの例を示す図。 コンピュータの構成を示す図。
　以下、本発明の実施の形態について、詳細に説明する。なお、同じ機能を有する構成部には同じ番号を付し、重複説明を省略する。
＜本発明とは異なる解決方法の候補１＞
　解決方法として、正常時の観測データと異常時の観測データを大量に用意して確率分布を生成しておく。図１は、正常時と異常時の確率分布の例である。対象の観測データが得られた際に、条件付き確率ｐ（ｘ｜ｙ＝１），ｐ（ｘ｜ｙ＝０）を求め、その確率密度比を
　Ｉｎ（ｐ（ｘ｜ｙ＝１）／ｐ（ｘ｜ｙ＝０））
のように求め、異常度とする方法があり得る。しかし、自動車内の状況は時々刻々と変化するため、正常パターン、異常パターンを網羅するほどの大量なデータを収集することは現実的ではない。したがって、正常・異常の確率分布を推定するのも困難である。よって、解決方法の候補１は本発明に利用できない。
＜本発明とは異なる解決方法の候補２＞
　図２に候補２の手法を示す。これは、非特許文献１に示された特異スペクトル変換を用いる手法である。過去と現在の確率分布の相違度（確率密度比）を、部分空間法を用いて評価する。特異スペクトル変換では、図２に示すように、過去および現在区間の時系列データの特異値分解により特徴的なパターンを抽出し、過去の特徴パターンが張る部分空間と現在の特徴パターンが張る部分空間の重なりの大きさから、変化度を計算する。なお、現在区間とは、変化度を求める対象の区間を意味している。具体的には、以下のように変化度を求める。
　図３にデータ行列の作成手順を示す。Ｋ，Ｍ，Ｎを２以上の整数、ｎを１以上Ｎ以下の整数、Ｌを１以上の整数、ｋを１以上Ｋ以下の整数、ｉ，ｒを１以上の整数、ｔを１以上Ｍ＋Ｎ－１以下の整数、［］は数値が行方向に並んだ行ベクトルまたは列ベクトルが行方向に並んだ行列を示す記号、Ｔは転置を示す記号、ｘ（ｉ）を時系列データとする。
　時系列データｘ（１），…，ｘ（Ｋ（Ｍ＋Ｎ－１））をＭ＋Ｎ－１個ずつの区間に分けて時系列データを要素とする区間Ｄ１，Ｄ２，…，Ｄｋ，…，ＤＫを生成する。Ｄｋは、ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）），…，ｘ（ｋ（Ｍ＋Ｎ－１））を区間データとする区間である。なお、Ｄｋは、
　Ｄｋ＝［ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）），…，ｘ（ｋ（Ｍ＋Ｎ－１））］Ｔ
のように、Ｍ＋Ｎ－１個の要素を持つ列ベクトルとして表現してもよい。区間データＤｋは、窓長（部分時系列の長さ）Ｍ、シフト幅１の部分時系列に分割され、Ｎ個の部分時系列が作成される。ｋ個目の区間Ｄｋのｎ番目の部分時系列ｘｋ，ｎは、Ｍ個の時系列データｘ（ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）），…，ｘ（ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）を要素とし、
　ｘｋ，ｎ＝［ｘ（ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）），…，ｘ（ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）］Ｔ
のように、列ベクトルで表現される。なお、表記に制限があるため、上記のように列ベクトルを表記しているが、この表記は、
Figure JPOXMLDOC01-appb-M000001
と同じ意味である。以下でも同じである。次に、区間ＤｋのＮ個の部分時系列ｘｋ，１，…，ｘｋ，Ｎをまとめて、データ行列Ｘｋを、
　Ｘｋ＝［ｘｋ，１…ｘｋ，ｎ…ｘｋ，Ｎ］
のように生成する。なお、表記に制限があるため、上記のように行列を表記しているが、この表記は、
Figure JPOXMLDOC01-appb-M000002
と同じである。ｘｋ，１，ｘｋ，ｎ，ｘｋ，Ｎの上下に延びた線は、ｘｋ，１，ｘｋ，ｎ，ｘｋ，Ｎが列ベクトルであることを示している。
　ここで、ｋ番目の区間Ｄｋを現在の区間（第二区間）とし、Ｌ個前の区間Ｄｋ－Ｌを比較対象の過去の区間（第一区間）とする。候補２での変化度の算出手法を図４に示す。なお、Ｌはラグと呼ばれ、過去区間（第一区間）と現在区間（第二区間）の相互位置を表す。第一区間および第二区間のデータ行列Ｘｋ－Ｌ，Ｘｋについて特異値分解を行い、特異値上位ｒ個の左特異ベクトルを取り出す。第一区間のデータ行列Ｘｋ－Ｌから抽出した特異ベクトルｕｋ－Ｌ，１，ｕｋ－Ｌ，２，…，ｕｋ－Ｌ，ｒの行列を以下で定義する。
　Ｕｋ－Ｌ，ｒ＝［ｕｋ－Ｌ，１，ｕｋ－Ｌ，２，…，ｕｋ－Ｌ，ｒ］
ただし、ｕｋ－Ｌ，１，ｕｋ－Ｌ，２，…，ｕｋ－Ｌ，ｒは列ベクトルである。同様に、第二区間のデータ行列Ｘｋ－Ｌから抽出した特異ベクトルｑｋ，１，ｑｋ，２，…，ｑｋ，ｒの行列を以下で定義する。
　Ｑｋ，ｒ＝［ｑｋ，１，ｑｋ，２，…，ｑｋ，ｒ］
ただし、ｑｋ，１，ｑｋ，２，…，ｑｋ，ｒは列ベクトルである。なお、ｒは第一区間と第二区間で異なる値を設定してもよい。
　変化度ｃ（ｋ）は、
　ｃ（ｋ）＝１－∥Ｕｋ－Ｌ，ｒ ＴＱｋ，ｒ∥２
ただし、∥　∥２は行列２ノルムを表し、∥Ｕｋ－Ｌ，ｒ ＴＱｋ，ｒ∥２はＵｋ－Ｌ，ｒ ＴＱｋ，ｒの最大特異値に等しい。
　実際の自動車内で観測された時系列データに対して特異値スペクトル変換を用いた変化検知結果を図５に示す。図５（Ａ）の縦軸は入力された時系列データの振幅を、図５（Ｂ）の縦軸は変化度を示しており、どちらも横軸は時間を示している。実験に用いた時系列データのサンプリング周波数は１６ｋＨｚであり、図５（Ａ）の矢印で示している範囲が、異常を検知すべき区間である。この実験では、部分時系列の長さＭ＝２５６、データ行列の列数Ｎ＝１２８、ラグＬ＝１２８とした。また、特異ベクトルの抽出数ｒは、特異値の累積寄与率が８０％を超えるまでの上位個数と設定した。図５（Ｂ）では、変化度のしきい値を０．４に設定し、変化度がしきい値を超えた点に三角の印を付し、異常と判断した点を強調している。図５より、異常を検出すべき区間以外にも異常と判断した点が散見されることが分かる。また、図５（Ａ）を確認すると、４３秒、８１秒、９９秒のあたりに大きな振動がある。図５（Ｂ）の結果では、しきい値を変更しても、この３か所だけを検出できるようにはできないことも分かる。
　非特許文献１の技術では、着目する区間における特異値分解により特徴パターンを抽出する。これは着目する区間のデータ行列の分散行列の固有値を分解していることに相当する。すなわち、着目する区間のみに着目し、分散の大きくなるような方向ベクトルを抽出しているため、着目区間に定常性を仮定したような、静的な特徴パターンしか抽出できていない。しかし、状況が時々刻々変化する環境では定常性を仮定できないため、特徴をうまく抽出できず、変化度の誤差が大きくなり、誤検出が多くなっていると考えられる。
＜本発明の異常推定装置＞
　図６に異常推定装置の機能構成例、図７に異常推定装置の処理フローの例を示す。図８に本発明の手法を示す。異常推定装置１０は、複数の利用者が運転する自動車において異常が発生した時刻を推定するための装置である。異常推定装置１０は、センサ２１０、第一取得部１１０、予測部１２０、異常度取得部１３０、異常判定部２２０を備える。センサ２１０は、自動車に設置され、振動に関する時系列データｘ（ｉ）を取得する（Ｓ２１０）。なお、以下では、第一区間（過去区間）は異常度を求めるための基準となる区間であり、第二区間（現在区間）は異常度を求める対象の区間である。
　第一取得部１１０は、時系列データｘ（ｉ）の一部区間である第一区間（過去区間）Ｄｋ－Ｌの動的な特徴パターンを取得する（Ｓ１１０）。動的な特徴パターンとは、周期振動と当該周期振動の減衰率に関する特徴パターンである。つまり、どのような周期振動で時間変化するか、その周期振動はどのような減衰（増加）率かを抽出する。より具体的には、動的モード分解によって、第一区間の動的な固有モードΦ（Ｍ次元の列ベクトルがランク数分並んだ行列）と固有値Λ＝ｄｉａｇ（λｊ）を求める。次に、動的モード分解について説明する。なお、正負の向きが変わるだけなので、本明細書においては、「減衰率」は増加する場合も含んでいる。
　図３を用いた候補２の説明と同じように、時系列データｘ（ｉ）をＭ＋Ｎ－１個ずつの区間に分けて時系列データを要素とする区間Ｄ１，Ｄ２，…，Ｄｋ，…，ＤＫとする。この中のｋ－Ｌ個目の区間Ｄｋ－Ｌが第一区間（過去区間）であり、ｋ個目の区間Ｄｋが第二区間（現在区間）である。第一取得部１１０は、区間Ｄｋ－Ｌについてのｎ番目の部分時系列ｘｋ－Ｌ，ｎが、
　ｘｋ－Ｌ，ｎ＝［ｘ（ｎ＋（ｋ－Ｌ－１）（Ｍ＋Ｎ－１）），…，ｘ（ｎ＋（ｋ－Ｌ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）］Ｔ
となるのように、列ベクトルである部分時系列ｘｋ，１，…，ｘｋ，Ｎを生成する。第一取得部１１０は、行列Ｘｋ－Ｌ，１＝［ｘｋ，１…ｘｋ，ｎ…ｘｋ，Ｎ－１］と行列Ｘｋ－Ｌ，２＝［ｘｋ，２…ｘｋ，ｎ…ｘｋ，Ｎ］を生成する。なお、行列Ｘｋ－Ｌ，１と行列Ｘｋ－Ｌ，２は、表記の関係で上記のように表現しているが、これらの表記は、
Figure JPOXMLDOC01-appb-M000003
Figure JPOXMLDOC01-appb-M000004
と同じ意味である。
　第一取得部１１０は、Ｘｋ－Ｌ，１を特異値分解した結果に基づいて求めた低ランク近似した疑似逆行列Ｘｋ－Ｌ，１ ＋を用いて、
　Ａ＝Ｘｋ－Ｌ，２Ｘｋ－Ｌ，１ ＋
のように線形作用素Ａを求める。そして、第一取得部１１０は、
　Ａ＝Φｋ－ＬΛｋ－ＬΦｋ－Ｌ －１
となる固有モードΦｋ－Ｌと固有値Λｋ－Ｌを求める（Ｓ１１０）。
　なお、特異値分解から疑似逆関数を用いて線形作用素Ａを求める詳細な手順は以下のとおりである。まずＸｋ－Ｌ，１を特異値分解すると、
　Ｘｋ－Ｌ，１＝Ｕｋ－Ｌ，１Σｋ－Ｌ，１Ｖｋ－Ｌ，１ Ｔ
のように特異行列Σｋ－Ｌ，１とユニタリ行列Ｕｋ－Ｌ，１，Ｖｋ－Ｌ，１が得られる。これらをランク数ｒまでで近似すると、
　Ｘｋ－Ｌ，１，ｒ＝Ｕｋ－Ｌ，１，ｒΣｋ－Ｌ，１，ｒＶｋ－Ｌ，１，ｒ Ｔ
となり、低ランク近似した疑似逆行列Ｘｋ－Ｌ，１ ＋は、
　Ｘｋ－Ｌ，１ ＋＝Ｖｋ－Ｌ，１，ｒΣｋ－Ｌ，１，ｒ －１Ｕｋ－Ｌ，１，ｒ Ｔ
のように求めることができる。したがって、線形作用素Ａは、Ｘｋ－Ｌ，１の特異値分解と低ランク近似によって、
　Ａ＝Ｘｋ－Ｌ，２Ｖｋ－Ｌ，１，ｒΣｋ－Ｌ，１，ｒ －１Ｕｋ－Ｌ，１，ｒ Ｔ
のように求めることができる。このように低ランク近似することで計算コストを削減できる。
　また、固有モードΦｋ－Ｌと固有値Λｋ－Ｌを求める詳細な手順は以下のとおりである。ここでもコスト削減のため低ランク近似を利用する。まず、
　～Ａ＝Ｕｋ－Ｌ，１，ｒ ＴＡＵｋ－Ｌ，１，ｒ
　　＝Ｕｋ－Ｌ，１，ｒ ＴＸｋ－Ｌ，２Ｖｋ－Ｌ，１，ｒΣｋ－Ｌ，１，ｒ －１
のように行列～Ａを求め、固有値分解により、
　～ＡＷ＝ＷΛ
となる固有ベクトルの列が並んだ行列Ｗと固有値Λｋ－Ｌ＝ｄｉａｇ（λｊ）を求める。そして、固有モードΦｋ－Ｌを
　Φｋ－Ｌ＝Ｘｋ－Ｌ，２Ｖｋ－Ｌ，１，ｒΣｋ－Ｌ，１，ｒ －１Ｗ
のように求めればよい。このように求めた固有モードΦｋ－Ｌと固有値固有値Λｋ－Ｌは、
　ＡΦｋ－Ｌ＝Φｋ－ＬΛ
の関係を満たしている。
　予測部１２０は、時系列データの一部区間である第一区間Ｄｋ－Ｌのデータである第一区間データに基づいて、第一区間Ｄｋ－Ｌよりも後の時系列データの一部区間である第二区間Ｄｋのデータを予測し、予測第二区間データを得る（Ｓ１２０）。具体的には、第二区間のデータを、第二区間の最初の部分時系列ｘｋ，１と特徴パターン（固有モードΦｋ－Ｌと固有値固有値Λｋ－Ｌ）を用いて予測し、予測第二区間データ＾ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）），…，＾ｘ（Ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）を得る。上述のように、第一取得部１１０によって、固有モードΦｋ－Ｌと固有値固有値Λｋ－Ｌが取得されている。予測部１２０は、ｎ＝２，…，Ｎについて、＾ｘｋ，ｎを
　＾ｘｋ，ｎ＝Φｋ－ＬΛｋ－Ｌ ｎ－１Φｋ－Ｌ －１ｘｋ，１
のように求める。または、予測部１２０は、まず
　ｂ＝Φｋ－Ｌ －１ｘｋ，１
を求め、ｎ＝２，…，Ｎについて、
　＾ｘｋ，ｎ＝Φｋ－ＬΛｋ－Ｌ ｎ－１ｂ
のように＾ｘｋ，ｎを求めてもよい。
　そして、予測部１２０は、部分時系列ｘｋ，１と予測した部分時系列＾ｘｋ，２，…，＾ｘｋ，Ｎ－１それぞれの第１要素を予測第二区間データの１番目からＮ－１番目までのデータ＾ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）），…，＾ｘ（Ｎ－１＋（ｋ－１）（Ｍ＋Ｎ－１））とし、予測した部分系列＾ｘｋ，Ｎの全要素を予測第二区間データのＮ番目からＭ＋Ｎ－１番目までのデータ＾ｘ（Ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）），…，＾ｘ（Ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）とすればよい。なお、表記の制限があるため“＾”は“ｘ”の前に記載しているが、図８に示すように“＾”を“ｘ”の上部に記載しているものと同じ意味である。上述のように予測第二区間データを求めるので、過去区間（第一区間）の動的な特徴パターンで、現在区間（第二区間）の動的なパターンがどの程度表現できるかを評価できるデータを得ることができる。
　異常度取得部１３０は、予測第二区間データと、時系列データにおける第二区間の実際のデータである第二区間データと、の差に基づき第二区間における異常度を取得する（Ｓ１３０）。つまり、異常度取得部１３０は、予測第二区間データ＾ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）），…，＾ｘ（Ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）と、時系列データにおける第二区間の実際のデータｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）），…，ｘ（Ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）から、異常度ｃ（ｋ）を
　ｃ（ｋ）
　＝（（＾ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１））－ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）））２
＋（＾ｘ（２＋（ｋ－１）（Ｍ＋Ｎ－１））－ｘ（２＋（ｋ－１）（Ｍ＋Ｎ－１）））２＋，…，＋（＾ｘ（ｔ＋（ｋ－１）（Ｍ＋Ｎ－１））－ｘ（ｔ＋（ｋ－１）（Ｍ＋Ｎ－１）））２＋，…，＋（（＾ｘ（ｋ（Ｍ＋Ｎ－１））－ｘ（ｋ（Ｍ＋Ｎ－１）））２）１／２
のように取得すればよい。なお、上記の記述は、
Figure JPOXMLDOC01-appb-M000005
と同じ意味である。
　異常判定部２２０は、取得された異常度に基づき第二区間において異常が発生しているかを判定する（Ｓ２２０）。例えば、異常度のしきい値を設け、しきい値を超える差があるときに異常が発生していると判断すればよい。
　実際の自動車内で観測された時系列データに対して、本発明の異常推定装置１０の手順で異常を推定した結果を図９に示す。図９（Ａ）の縦軸は入力された時系列データの振幅を、図９（Ｂ）の縦軸は異常度を示しており、どちらも横軸は時間を示している。図９（Ａ）のデータは、図５（Ａ）と同じである。図５で示した実験と同じように、実験に用いた時系列データのサンプリング周波数は１６ｋＨｚであり、図９（Ａ）の矢印で示している範囲が、異常を検知すべき区間である。この実験では、部分時系列の長さＭ＝２５６、履歴行列およびテスト行列の列数Ｎ＝１２８、ラグＬ＝１２８とした。また、ランク数ｒは、特異値の累積寄与率が８０％を超えるまでの上位個数と設定した。図９（Ｂ）では、異常度のしきい値を０．３に設定し、異常度がしきい値を超えた点に三角の印を付し、異常と判断した点を強調している。このしきい値の設定が、上述の異常判定部２２０に相当する。図９（Ｂ）では、異常を検出すべき区間に集中的に異常と判断した点あることが分かる。特に、図９（Ｂ）では４３秒、８１秒、９９秒のあたりにある大きめな振動だけを検出できている。したがって、本発明の異常推定装置によれば、異常がない範囲を除外できるので、複数の利用者が運転する自動車に付いた傷の原因を調査する際に、確認しなければならない映像を限定できる。
＜本発明の異常検出装置＞
　実施例１は、シェアカー、レンタカーなどの複数の利用者が運転する自動車に用途を限定していた。しかし、本発明は、時間の経過に応じて特徴パターンが変化する可能性のある時系列データに対して幅広く適用できる。つまり、センサが取得するデータは振動に関するデータに限る必要はない。そこで、実施例２では、何らかの時系列データを入力とし、異常検出しやすい異常度を出力する異常検出装置について説明する。
　図１０に本発明の異常検出装置の機能構成例、図１１に本発明の異常検出装置の処理フローの例を示す。異常検出装置２０は、時系列データにおける異常度を検出するための装置である。異常検出装置２０は、第一取得部１１０、予測部１２０、異常度取得部１３０を備える。実施例１の異常推定装置１０との関係を示すため、図１０では異常推定装置１０、センサ２１０、異常判定部２２０を点線で示している。
　異常検出装置２０には、時系列データｘ（ｉ）が入力される（Ｓ２１１）。第一取得部１１０は、時系列データの一部区間である第一区間（過去区間）Ｄｋ－Ｌから第一区間の動的な特徴パターンを取得する（Ｓ１１０）。予測部１２０は、時系列データの一部区間である第一区間Ｄｋ－Ｌのデータである第一区間データに基づいて、第一区間Ｄｋ－Ｌよりも後の時系列データの一部区間である第二区間Ｄｋのデータを予測し、予測第二区間データを得る（Ｓ１２０）。異常度取得部１３０は、予測第二区間データと、時系列データにおける第二区間の実際のデータである第二区間データと、の差に基づき第二区間における異常度を取得する（Ｓ１３０）。ステップＳ１１０，Ｓ１２０，Ｓ１３０は、実施例１と同じである。
　実施例１で候補２として説明したような過去区間の固有モードと現在区間の固有モードの重なり度合いから変化度を求める手法では、減衰（増加）率を表す固有値が考慮されていないため、同一の振動数を持った固有モードがあったとしても、過去区間では減衰する特徴パターン、現在区間では増加する特徴パターンに相当していた場合、固有モードだけの比較では区別がつかないため、異常検出に利用できるような変化度を取得できなかった。
　本発明によれば、時系列データのみから、注目する区間の動的な特徴パターンを推定できるため、状況が動的に変化する系から得られる時系列データの変化度を推定する精度が向上し、変化部位の検出精度が向上する。観測区間で抽出された特徴パターンは、固有モード（固有振動数）に対応した信号であるため、データの解釈が容易になる。したがって、本発明の異常検出装置によれば、動的な特徴パターンに基づいて予測第二区間データを得るので、時間の経過に応じて時系列データの特徴パターンが変化したとしても、異常検出に必要な異常度を取得できる。
［プログラム、記録媒体］
　上述の各種の処理は、図１２に示すコンピュータの記録部２０２０に、上記方法の各ステップを実行させるプログラムを読み込ませ、制御部２０１０、入力部２０３０、出力部２０４０などに動作させることで実施できる。
　この処理内容を記述したプログラムは、コンピュータで読み取り可能な記録媒体に記録しておくことができる。コンピュータで読み取り可能な記録媒体としては、例えば、磁気記録装置、光ディスク、光磁気記録媒体、半導体メモリ等どのようなものでもよい。
　また、このプログラムの流通は、例えば、そのプログラムを記録したＤＶＤ、ＣＤ－ＲＯＭ等の可搬型記録媒体を販売、譲渡、貸与等することによって行う。さらに、このプログラムをサーバコンピュータの記憶装置に格納しておき、ネットワークを介して、サーバコンピュータから他のコンピュータにそのプログラムを転送することにより、このプログラムを流通させる構成としてもよい。
　このようなプログラムを実行するコンピュータは、例えば、まず、可搬型記録媒体に記録されたプログラムもしくはサーバコンピュータから転送されたプログラムを、一旦、自己の記憶装置に格納する。そして、処理の実行時、このコンピュータは、自己の記録媒体に格納されたプログラムを読み取り、読み取ったプログラムに従った処理を実行する。また、このプログラムの別の実行形態として、コンピュータが可搬型記録媒体から直接プログラムを読み取り、そのプログラムに従った処理を実行することとしてもよく、さらに、このコンピュータにサーバコンピュータからプログラムが転送されるたびに、逐次、受け取ったプログラムに従った処理を実行することとしてもよい。また、サーバコンピュータから、このコンピュータへのプログラムの転送は行わず、その実行指示と結果取得のみによって処理機能を実現する、いわゆるＡＳＰ（Application Service Provider）型のサービスによって、上述の処理を実行する構成としてもよい。なお、本形態におけるプログラムには、電子計算機による処理の用に供する情報であってプログラムに準ずるもの（コンピュータに対する直接の指令ではないがコンピュータの処理を規定する性質を有するデータ等）を含むものとする。
　また、この形態では、コンピュータ上で所定のプログラムを実行させることにより、本装置を構成することとしたが、これらの処理内容の少なくとも一部をハードウェア的に実現することとしてもよい。
１０　異常推定装置
２０　異常検出装置
１１０　第一取得部
１２０　予測部
１３０　異常度取得部
２１０　センサ
２２０　異常判定部

Claims (8)

  1. 　複数の利用者が運転する自動車において異常が発生した時刻を推定するための異常推定装置であって、
    　前記自動車に設置され、振動に関する時系列データを取得するセンサと、
    　前記時系列データの一部区間である第一区間のデータである第一区間データに基づいて、前記第一区間よりも後の前記時系列データの一部区間である第二区間のデータを予測し、予測第二区間データを得る予測部と、
    　前記予測第二区間データと、前記時系列データにおける第二区間の実際のデータである第二区間データと、の差に基づき前記第二区間における異常度を取得する異常度取得部と、
    　取得された異常度に基づき前記第二区間において異常が発生しているかを判定する異常判定部と、
    を有する異常推定装置。
  2. 　請求項１記載の異常推定装置であって、
    　前記第一区間の動的な特徴パターンを取得する第一取得部も有し、
    　前記予測部は、前記特徴パターンを用いて前記予測第二区間データを予測する
    　ことを特徴とする異常推定装置。
  3. 　時系列データにおける異常度を検出する異常検出装置であって、
    　前記時系列データの一部区間である第一区間から前記第一区間の動的な特徴パターンを取得する第一取得部と、
    　前記第一区間よりも後の前記時系列データの一部区間である第二区間のデータを、前記特徴パターンを用いて予測し、予測第二区間データを得る予測部と、
    　前記予測第二区間データと、前記時系列データにおける前記第二区間の実際のデータと、の差に基づき異常度を取得する異常度取得部と、
    を有する異常検出装置。
  4. 　時系列データにおける異常度を検出する異常検出装置であって、
    　前記時系列データの一部区間である第一区間における周期振動と当該周期振動の減衰率に関する特徴パターンを取得する第一取得部と、
    　前記第一区間よりも後の前記時系列データの一部区間である第二区間のデータを、前記特徴パターンを用いて予測し、予測第二区間データを得る予測部と、
    　前記予測第二区間データと、前記時系列データにおける前記第二区間の実際のデータと、の差に基づき異常度を取得する異常度取得部と、
    を有する異常検出装置。
  5. 　時系列データにおける異常度を検出する異常検出装置であって、
    　前記時系列データの一部区間である第一区間のデータの動的モード分解によって、固有モードと固有値を求める第一取得部と、
    　前記第一区間よりも後の前記時系列データの一部区間である第二区間のデータを、当該第二区間の最初のあらかじめ定めた数のデータと前記固有モードと前記固有値を用いて予測し、予測第二区間データを得る予測部と、
    　前記予測第二区間データと、前記時系列データにおける前記第二区間の実際のデータと、の差に基づき異常度を取得する異常度取得部と、
    を有する異常検出装置。
  6. 　時系列データにおける異常度を検出する異常検出装置であって、
    　ＭとＮを２以上の整数、ｎを１以上Ｎ以下の整数、Ｌを１以上の整数、ｋをＬ＋１以上の整数、ｉを１以上の整数、ｔを１以上Ｍ＋Ｎ－１以下の整数、［］は数値が行方向に並んだ行ベクトルまたは列ベクトルが行方向に並んだ行列を示す記号、Ｔは転置を示す記号、ｘ（ｉ）を時系列データ、ｘｋ，ｎをｘｋ，ｎ＝［ｘ（ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）），…，ｘ（ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）］Ｔである列ベクトル、Ｘｋ－Ｌ，１をＸｋ－Ｌ，１＝［ｘｋ－Ｌ，１…ｘｋ－Ｌ，ｎ…ｘｋ－Ｌ，Ｎ－１］である行列、Ｘｋ－Ｌ，２をＸｋ－Ｌ，２＝［ｘｋ－Ｌ，２…ｘｋ－Ｌ，ｎ…ｘｋ－Ｌ，Ｎ］である行列とし、
    　前記時系列データの一部区間である第一区間のデータｘ（１＋（ｋ－Ｌ－１）（Ｍ＋Ｎ－１）），…，ｘ（Ｎ＋（ｋ－Ｌ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）において、
    　Ｘｋ－Ｌ，１を特異値分解した結果に基づいて求めた低ランク近似した疑似逆行列Ｘｋ－Ｌ，１ ＋を用いて、
    　Ａ＝Ｘｋ－Ｌ，２Ｘｋ－Ｌ，１ ＋
    のように線形作用素Ａを求め、
    　Ａ＝Φｋ－ＬΛｋ－ＬΦｋ－Ｌ －１
    となる固有モードΦｋ－Ｌと固有値Λｋ－Ｌを求める第一取得部と、
    　前記第一区間よりも後の前記時系列データの一部区間である第二区間のデータを、
    　＾ｘｋ，ｎ＝Φｋ－ＬΛｋ－Ｌ ｎ－１Φｋ－Ｌｘｋ，１
    　のように予測し、予測第二区間データ＾ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）），…，＾ｘ（Ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）を得る予測部と、
    　前記予測第二区間データ＾ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）），…，＾ｘ（Ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）と、前記時系列データにおける前記第二区間の実際のデータｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）），…，ｘ（Ｎ＋（ｋ－１）（Ｍ＋Ｎ－１）＋Ｍ－１）から、異常度ｃ（ｋ）を
    　ｃ（ｋ）
    　＝（（＾ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１））－ｘ（１＋（ｋ－１）（Ｍ＋Ｎ－１）））２
    ＋（＾ｘ（２＋（ｋ－１）（Ｍ＋Ｎ－１））－ｘ（２＋（ｋ－１）（Ｍ＋Ｎ－１）））２＋，…，＋（＾ｘ（ｔ＋（ｋ－１）（Ｍ＋Ｎ－１））－ｘ（ｔ＋（ｋ－１）（Ｍ＋Ｎ－１）））２＋，…，＋（（＾ｘ（ｋ（Ｍ＋Ｎ－１））－ｘ（ｋ（Ｍ＋Ｎ－１）））２）１／２
    のように取得する異常度取得部と、
    を有する異常検出装置。
  7. 　第一取得部と予測部と異常度取得部を備えた異常検出装置を用いて、時系列データにおける異常度を検出する異常検出方法であって、
    　前記第一取得部が、前記時系列データの一部区間である第一区間から前記第一区間の動的な特徴パターンを取得する第一取得ステップと、
    　前記予測部が、前記第一区間よりも後の前記時系列データの一部区間である第二区間のデータを、前記特徴パターンを用いて予測し、予測第二区間データを得る予測ステップと、
    　前記異常度取得部が、前記予測第二区間データと、前記時系列データにおける前記第二区間の実際のデータと、の差に基づき異常度を取得する異常度取得ステップと、
    を実行する異常検出方法。
  8. 　請求項３から６のいずれかに記載の異常検出装置としてコンピュータを機能させるための異常検出プログラム
PCT/JP2019/030717 2019-08-05 2019-08-05 異常推定装置、異常検出装置、異常検出方法、異常検出プログラム WO2021024351A1 (ja)

Priority Applications (3)

Application Number Priority Date Filing Date Title
JP2021538566A JP7248127B2 (ja) 2019-08-05 2019-08-05 異常推定装置、異常検出装置、異常検出方法、異常検出プログラム
PCT/JP2019/030717 WO2021024351A1 (ja) 2019-08-05 2019-08-05 異常推定装置、異常検出装置、異常検出方法、異常検出プログラム
US17/632,527 US12090940B2 (en) 2019-08-05 2019-08-05 Abnormality estimating apparatus, abnormality detection apparatus, abnormality detection method, and abnormality detection program

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
PCT/JP2019/030717 WO2021024351A1 (ja) 2019-08-05 2019-08-05 異常推定装置、異常検出装置、異常検出方法、異常検出プログラム

Publications (1)

Publication Number Publication Date
WO2021024351A1 true WO2021024351A1 (ja) 2021-02-11

Family

ID=74502891

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2019/030717 WO2021024351A1 (ja) 2019-08-05 2019-08-05 異常推定装置、異常検出装置、異常検出方法、異常検出プログラム

Country Status (3)

Country Link
US (1) US12090940B2 (ja)
JP (1) JP7248127B2 (ja)
WO (1) WO2021024351A1 (ja)

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH08320721A (ja) * 1995-05-24 1996-12-03 Meidensha Corp 故障予知装置の異常内容判定方法
JP2005094361A (ja) * 2003-09-17 2005-04-07 Intelligent Cosmos Research Institute ネットワーク異常検出装置、ネットワーク異常検出方法およびネットワーク異常検出プログラム
JP2005128808A (ja) * 2003-10-24 2005-05-19 Mitsubishi Electric Corp 予測装置
JP2006135412A (ja) * 2004-11-02 2006-05-25 Tokyo Gas Co Ltd 遠隔監視システム
JP2009111537A (ja) * 2007-10-29 2009-05-21 Nippon Telegr & Teleph Corp <Ntt> 通信関係構造変化による異常検出方法、装置、プログラム、および記録媒体

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US7177738B2 (en) 2001-05-30 2007-02-13 Alpine Electronics, Inc. Vehicle management system
JP4371616B2 (ja) * 2001-09-18 2009-11-25 本田技研工業株式会社 車両用衝突検出装置
GB201407952D0 (en) * 2014-01-31 2014-06-18 Cambridge Consultants Monitoring device
GB201301710D0 (en) 2013-01-31 2013-03-20 Cambridge Consultants Condition Monitoring Device

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JPH08320721A (ja) * 1995-05-24 1996-12-03 Meidensha Corp 故障予知装置の異常内容判定方法
JP2005094361A (ja) * 2003-09-17 2005-04-07 Intelligent Cosmos Research Institute ネットワーク異常検出装置、ネットワーク異常検出方法およびネットワーク異常検出プログラム
JP2005128808A (ja) * 2003-10-24 2005-05-19 Mitsubishi Electric Corp 予測装置
JP2006135412A (ja) * 2004-11-02 2006-05-25 Tokyo Gas Co Ltd 遠隔監視システム
JP2009111537A (ja) * 2007-10-29 2009-05-21 Nippon Telegr & Teleph Corp <Ntt> 通信関係構造変化による異常検出方法、装置、プログラム、および記録媒体

Also Published As

Publication number Publication date
JP7248127B2 (ja) 2023-03-29
JPWO2021024351A1 (ja) 2021-02-11
US20220289131A1 (en) 2022-09-15
US12090940B2 (en) 2024-09-17

Similar Documents

Publication Publication Date Title
US10914608B2 (en) Data analytic engine towards the self-management of complex physical systems
US10496084B2 (en) Dequantizing low-resolution IoT signals to produce high-accuracy prognostic indicators
US20180308012A1 (en) Method and device for improving the robustness against &#34;adversarial examples&#34;
US11573877B2 (en) Distributed architecture for fault monitoring
Dervilis et al. Robust methods of inclusive outlier analysis for structural health monitoring
US8347144B2 (en) False alarm mitigation
JP6402541B2 (ja) 異常診断装置及びプログラム
JP2004531815A (ja) 予測的状態監視のための診断システムおよび方法
JP2015170121A (ja) 異常診断装置及びプログラム
Ma et al. Deep recurrent convolutional neural network for remaining useful life prediction
Gres et al. Orthogonal projection-based harmonic signal removal for operational modal analysis
US11093863B2 (en) Method for ascertaining a time characteristic of a measured variable, prediction system, actuator control system, method for training the actuator control system, training system, computer program, and machine-readable storage medium
Stull et al. On assessing the robustness of structural health monitoring technologies
WO2021019672A1 (ja) 異常度推定装置、異常度推定方法、およびプログラム
CN117063173A (zh) Ai模型水印系统和方法
Lakshmi et al. Structural damage detection using ARMAX time series models and cepstral distances
WO2021024351A1 (ja) 異常推定装置、異常検出装置、異常検出方法、異常検出プログラム
KR102307076B1 (ko) 딥 앙상블과 칼만필터를 통합한 방식의 차량 사이드슬립 앵글 추정 방법 및 장치
JP7058386B2 (ja) 分析装置、分析方法およびプログラム
Mendler et al. On the probability of localizing damages based on mode shape changes
WO2020158032A1 (ja) 推定装置、学習装置、それらの方法、およびプログラム
Afsharmovahed et al. A Novel Damage Detection Approach Based on Feature Extraction and Selection Using Machine Learning Without Signal Processing: A Case Study on the Tianjin Yonghe Bridge
Lagerblad A study of vehicle fatigue damage predictions based on state estimations
US12073250B2 (en) Determining memory requirements for large-scale ml applications to facilitate execution in GPU-embedded cloud containers
US20230035541A1 (en) Optimizing a prognostic-surveillance system to achieve a user-selectable functional objective

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

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2021538566

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

Country of ref document: EP

Kind code of ref document: A1