WO2018122585A1 - Method for urban road traffic incident detecting based on floating-car data - Google Patents

Abstract

A method for urban road traffic incident detecting based on floating-car data can use on-board GNSS positioning devices of floating cars to acquire the space position information of the floating cars at different time, and therefore, can achieve the intelligent detection of urban road traffic incidents by analyzing and mining the massive floating car track information. The urban road traffic incident detection technology uses the probability distribution of travel speed to characterize traffic conditions, uses the measurement index of difference in probability distribution to reflect the difference in traffic conditions, thereby having the advantages of clear principles, being easy to implement, and having a high detection rate.

Description

Urban road traffic anomaly detection method based on floating car data

 A Method for Urban Traffic Incident Detecting based on Floating-Car Data

 The invention belongs to the technical field of traffic detection. In particular, the present invention relates to a method for real-time detection of urban road traffic anomalies. Through the on-board GNSS positioning device of the floating car, the spatial position information of different time points can be obtained. After data preprocessing, map matching and data fusion, the travel speed probability distribution of the specific time and space range can be obtained; according to the change of the speed distribution, the position can be effectively identified. Urban road traffic anomalies. Background technique

 Traffic anomaly detection is an important part of urban traffic management and one of the core functions of intelligent transportation systems. Traffic anomalies mainly include traffic accidents, vehicle breakdown, truck falling, damage or malfunction of road traffic facilities, and other special events that cause traffic flow disturbances. Such incidents are prone to traffic congestion, reduced road capacity, and severely affect the normal operation of the entire road traffic system. Through traffic anomaly detection, traffic managers can timely understand traffic anomaly information and take appropriate inducement and control measures to reduce the adverse effects of traffic anomalies.

 Traffic anomaly detection can be divided into manual mode and automatic mode. Manual methods include patrol cars, emergency telephone reporting and video surveillance. Due to the consumption of manpower and material resources and poor real-time performance, traffic management needs cannot be met. The automatic method relies on the automatic event detection (AID, Automated Incidence Detection) algorithm. The basic principle is to identify traffic anomalies by detecting changes in road traffic at different locations. Currently used AID algorithms include pattern recognition algorithms (such as California algorithm, Monica algorithm), statistical prediction algorithms (such as exponential smoothing, Kalman filtering), traffic flow model algorithms (such as McMaster algorithm), and intelligent recognition algorithms. (such as artificial neural networks, fuzzy logic algorithms).

 However, the current detection methods have the disadvantages of high requirements on facilities, high computational complexity, and inability to make further judgments on the situation of abnormal conditions. The invention utilizes the trajectory data returned by the taxi and the bus GNSS positioning device, establishes a historical traffic state database and a real-time traffic state database, and identifies the traffic anomaly event by analyzing the difference of the traffic flow characteristics reflected by the two. The method has the characteristics of good real-time performance, parallel processing, high recognition rate and low requirements for detection facilities. It is suitable for the detection of urban road traffic anomalies in data environment with real-time floating vehicle positioning data.

 At present, for the monitoring of traffic anomalies, the following representative technologies are available:

 A US patent application, US 20160148512, discloses a composition principle and implementation method of a traffic anomaly detection and reporting system. The system consists of a sensor, a communication module, a mobile processing module, and a user interaction module. The sensor is used to collect relevant data around the vehicle; the communication module is used for transmitting the vehicle data and receiving data of the surrounding vehicle; the mobile processing module is for processing and analyzing the data of the relevant vehicle in a certain area and generating a traffic event report; user interaction The module is able to provide traffic incident reports like a user. The scheme is a traffic anomaly detection technology based on the vehicle and vehicle communication network, which can use various types of information collected by the sensor to identify abnormal events. However, since the sensor and the communication unit need to be separately installed and debugged, the implementation is difficult; the processing capacity of the mobile processing unit is limited; and the mobile and fixed message receiving end is required, and the system itself has a failure probability and the reliability is not good.

A Chinese patent application, CN 104809878 A, discloses a method for detecting abnormal state of urban road traffic using bus GPS data. The scheme obtains the link delay time index according to the GPS historical data, obtains the instantaneous speed, the cycle average speed, the weighted moving average speed and the multi-vehicle average speed according to the current GPS data, and uses the gauge variable analysis algorithm to detect the abnormality. This program does not New testing facilities are needed and implementation is convenient. However, the characterization of the traffic situation is too simplistic, and it is impossible to analyze the characteristics and causes of traffic anomalies. There is no basis for the division of traffic scenarios, and the influence of weather and other factors on traffic situation changes cannot be considered. Summary of the invention

 In order to explain the contents of the present invention more clearly, the technical terms involved will first be explained as follows:

 Floating car: Also known as the probe car. Refers to buses and taxis that have on-board positioning devices and are driving on city roads.

 GNSS: Global Navigation Satellite System „ includes GPS, Glonass, Galileo and Beidou satellite navigation systems.

 Space-time sub-zone: A zone divided by two dimensions, time and space, reflected in a certain space within a certain period of time.

 Historical trajectory data: Historical trajectory data is trajectory data accumulated over a long period of time and stored in a database. Historical trajectory data is dynamic change data that needs to be updated in a timely manner and periodically reprocessed and analyzed to ensure the accuracy of historical traffic feature extraction. The data for each spatiotemporal sub-area can be processed in parallel to increase efficiency. In the present invention, it may be simply referred to as historical data.

 Real-time trajectory data: The real-time trajectory data is a trajectory data set within a time zone that is closest to the current time. In the present invention, it may be simply referred to as real-time data.

 Traffic situation: A general term for the comprehensive situation of traffic operations within a certain period of time and within a certain space.

 Traffic anomalies: traffic turbulence caused by events such as traffic accidents, vehicle breakdown, truck falling, road traffic facilities damage or malfunctions.

 Abnormal traffic severity: The severity of traffic flow disorder is the difference in traffic flow characteristics after traffic flow and traffic anomalies in normal conditions.

 Traffic Anomaly Index: A measure of the severity of traffic. The range is 0~10. The larger the value, the more serious the traffic anomaly.

 Traffic environment: The sum of all external influences and forces acting on road traffic participants. This includes road conditions, transportation facilities, landforms, meteorological conditions, and transportation activities of other transportation participants.

 Map Matching: The process of associating geographic coordinates with a city road network.

 Peak hourly traffic: The maximum hourly traffic flow in a city's road section.

Finite Mixing Model: A mathematical method of simulating complex density with simple density. A finite mixed model with a set of variables y and a component number K can be expressed as:

 Response variable: A variable that changes according to an independent variable, also called a dependent variable.

 Bayesian information criterion: It is an evaluation index of the reliability of the result of correcting the probability of occurrence by using the Bayesian formula to estimate the probability of partial unknowns under incomplete intelligence. Its calculation method is:

 Where, is the maximum value of the likelihood function, which is the number of unknown parameters, "for the sample size.

 Likelihood function: A likelihood function is a function of the parameters of a statistical model. Given the output X, the likelihood function L(0\x) on the parameter Θ (in numerical value) is equal to the probability of the variable after the given parameter: ( |χ)=Ρ( =χ| ).

 Kullback-Leibler divergence: A measure of the difference between two probability distributions P and Q.

Jensen-Shannon divergence: is a symmetrized form of Kullback-Leibler divergence. K-Medoids algorithm: A clustering algorithm that selects such a point from the current category for each iteration to all other

The sum of the distances of points (in the current category) is the smallest as a new center point. The object of the present invention is to establish a scheme based on a floating vehicle trajectory recording system, using historical GNSS positioning data and real-time GNSS positioning data, combined with traffic environment information to identify road traffic anomalies. In order to achieve the above object, the present invention provides the following technical solutions:

 The implementation premise of the present invention is: a floating car (a taxi, a bus, etc.) equipped with a GNSS track recorder; a data center having a large-scale storage, calculation, and real-time task processing capability.

 The scope of application of the present invention is: Urban roads (including ground roads and elevated roads) through which the above-mentioned floating vehicles pass.

 The implementation steps of the present invention include:

 1) Establish a space-time sub-zone. The day is divided into time segments, and the implementation area of urban road traffic anomaly detection is divided into several spatial segments.

 2) Data preprocessing. Data cleaning, data integration, data conversion, and data reduction of GNSS positioning data to improve the structure of data.

 3) Quick map matching. Combined with the urban road network data, the map matching algorithm is used to project the GNSS positioning points to the city map, and the matching relationship between the positioning points and the road segments is established, and the error caused by the positioning drift is corrected.

 4) Data sampling. Different sampling methods can be used based on different computing power and accuracy requirements.

 5) Historical trajectory data analysis and feature extraction. Using the historical floating vehicle trajectory data, a traffic feature model is established. For each spatiotemporal sub-area, the traffic characteristics are described by the probability distribution of the travel speed, the traffic feature model is established and the parameters are estimated.

 6) Real-time trajectory data analysis and feature extraction. Use real-time floating car trajectory data to grasp the dynamics of traffic characteristics. The current traffic characteristics are described by the travel speed probability distribution of the current space-time sub-region, and the model is established and parameter estimation is performed.

 7) Anomaly detection. The Jensen-Shannon divergence is used to measure the difference between historical traffic characteristics and real-time traffic characteristics.

 8) Quantitative characterization of abnormal severity. Calculate and publish an abnormality index for traffic conditions. The step 1) may specifically adopt the following methods:

 11) Isometric space-time division method. The segment size of the time dimension is determined. The time segment span is a fixed value, usually 30 minutes as a time segment; the segment size of the spatial dimension is determined, and the spatial segment span is a fixed value, and a spatial grid of 200 m×200 m is usually taken as a spatial segment.

12) Non-equidistant space-time division method. For urban central areas where the road network density is greater than 2km/km ² or the peak hour flow is greater than 1000 vehicles/hour, 30min time segments and 200mX 200m space segments are taken. For road network density less than 2km/km ² or peak hour traffic is less than 1000. In the suburb of the city/hour, take a 30-minute time segment and a 400mX 400m space segment. The step 3) specifically includes the following steps:

 31) Divide the space area to be processed into a grid of a certain size, and the range of each grid area can be expressed as

32) Determine the grid area where the anchor point is located, and use the distance and azimuth to search for the road segment where the anchor point is located. The matching scheme includes:

321) Single point matching scheme: Searching for the nearest road segment from point A, when the difference between the traveling direction angle satisfying point A and the direction angle of the road segment ij is less than the threshold value, that is, |-.|< is satisfied, and the matching is completed; if not satisfied, - -.|< Delete the road segment zy in the search space and continue to search for other road segments until the conditions are met. The matching method is shown in Figure 3.

 322) Point sequence matching scheme:

This program is suitable for high frequency floating car data. The floating vehicle GNSS data acquisition frequency is expressed as /. =1/3⁄4, the point P( -t ₀ ), P( _A +t _Q ) adjacent to A in time is defined as 1-adjacent point of A, ■Pfe o), Pt _A +2k, defined as A 2-adjacent points, and so on, then corpse 04-), ·Ρ(ί4+) is defined as the neighboring point of A. in/. When <1Ηζ, take ^=1 or 2. Take the road segment zy with the smallest distance from the neighboring points of A and A, and calculate the mean value of the driving direction angle of the 1-adjacent point of A and A^41. If it satisfies |1|<, complete the matching; otherwise, search for other road segments until it is satisfied. condition.

 33) Use the straight line equation of the road segment (if it is a curved road segment, it is roughly split into straight lines), calculate the projection coordinates of the GNSS positioning point on the road segment, and reduce the error caused by the GNSS positioning drift. The specific method uses the GNSS positioning point linear projection method as:

Determine the straight line equation of the road segment zy (if the road segment is a curve, divide it into several straight line segments): yy _i =k(xx _i ) where the slope is:

XJ - ^χ , J The projection line equation is: yy _A = - ( ^X - ^X A) ky _A -ky. + k ² x. + x _A

Solve the projected coordinate p as: ^Xp= t^r" ^

k ² y _A + y _{ + hc _A - hc _i After the map matching process, combined with the timestamp data of the coordinates of the positioning point, the positioning point is matched to the space-time sub-region. The step 4) may specifically adopt the following methods:

41) Full sample program. The data of all driving speeds of each secondary floating car in a time and space sub-region constitutes the whole. The implementation method is to calculate the travel speed of each vehicle in the space-time sub-region: = ^{2 + 3} two... ""- ¹ '" , where ^...Α- is the first and second GNSS positioning in the space-time sub-region The distance between points, ..., the distance between the n-1th and the nth GNSS anchor point, h-t„ is the first in the space-time sub-region, ..., Timestamps of GNSS anchor points; data in each spatio-temporal sub-area is not filtered to form a set for subsequent processing.

42) Time-smooth sampling plan. Specify the length of the time segment and the upper limit of the number of segments of the same time; search for the velocity data in each time segment in a time-space sub-region. If the number of velocity data in the time segment exceeds the upper limit, the data of the upper limit is randomly added to the data to be processed. sample. The implementation method is to calculate the travel speed of each vehicle in the space-time sub-area: V- = ² + ³ +... + - , Where Α, ₂ ... ί . is the distance between the first and second GNSS anchor points in the space-time sub-region, ..., the first -1 and the first

The distance between the GNSS positioning points is the first time in the space-time sub-area, ..., the time stamp of the GNSS positioning point; the specified time segment length t _P , the upper limit of the number of segments of the same time; 3⁄4. _1; Search for the speed data in the time and time segments of a time-space sub-region, if the number of speed data in the time segment exceeds the upper limit; 1⁄2 _∞ , random take; 1⁄2 _∞ data is added. The step 5) may specifically adopt the following method:

 51) Simple historical trajectory data fusion method. The historical data under the condition of no traffic anomalies, as a whole, the traffic characteristic model establishment and parameter estimation. The method uses a finite mixing model to establish a traffic feature model and perform parameter estimation. One of the following three options can be used:

 511) Mixed Gaussian model of fixed composition

 This scheme uses a mixed Gaussian model with a fixed component quantity to describe the probability distribution of vehicle speed. The number of components is manually specified according to the distribution pattern of the vehicle speed under typical conditions. In order to ensure the reliability of the probability distribution, the number of components f cannot be too small. Generally, f=4~6 is available.

 512) Mixed Gaussian model with variable composition

 This program uses a model-based evaluation method to select the appropriate number of components, as follows:

 Determine the maximum number of possible components K, and estimate the parameters of the mixed Gaussian model of “=1, 2,... components separately; for each model, determine the best model by Bayesian information criterion (.BIC,). The maximum number of components is generally selected according to the accuracy requirement, but it must be noted that the more components, the slower the expectation maximization algorithm converges. The maximum number of components selected here is ^=5, which requires calculation:

/ ^ν ) = ί Λ( ^ν ) " ε {ΐ, 2".., 5}

 /=1

 There are 5 mixed models in total. At the same time, the definition of the five models is defined as:

 BIC = -2\n L + k - \n n

 Where, is the maximum likelihood function value, t is the number of parameters in the model, "for the total amount of data.

 After that, the /C smallest hybrid model is selected, and its parameter vectors 1], μ, and σ are recorded as the feature records of the local space-time sub-region. The density curve morphology of the hybrid model is shown in Figure 6.

 513) Finite mixed model with variable composition and distribution type

 This scheme uses the same model-based evaluation method as 512), but the distribution of the sub-components and the number of components are variable, as follows:

 The probability distribution model is chosen as the distribution type of the sub-components, including but not limited to: normal distribution, gamma distribution, Weibull distribution. When using a normal distribution, the sub-distribution function takes:

 ,2,

 - y)

Feet: -exp When using a gamma distribution, the sub-distribution function takes: =^ ^ ^a "- ^le ' where ^Γ( ) = ^ ' ^dt

When using a Weibull distribution, the sub-distribution function takes:

 Assuming that the distribution types of all sub-components of the hybrid model are the same, determine the maximum number of possible components ^. For the selection of the M seed component distribution type and the number of f species, a total of M combinations are formed, and the /C value is calculated separately, and the model with the smallest /C is the best model.

 52) Historical trajectory data classification method by context. Based on temperature, precipitation, visibility and traffic control measures, historical data under no traffic anomalies are divided into different categories, models are established and parameter estimates are made. The implementation method is as follows:

 According to the temperature, precipitation, visibility and traffic control measures, the traffic environment is divided into 5~8 categories, and the historical data is classified into the above categories according to the different traffic environments corresponding to historical data. For each category, the processing as described in 51) is performed separately, thereby establishing a mapping relationship 3⁄4 →7, which is a traffic situation and a traffic situation.

 53) Historical data clustering method. For the historical data, the difference quantization between different spatiotemporal regions is obtained by comparison between the time and space sub-regions, and the quantized differences are used for clustering. Temperature, precipitation, visibility and traffic control measures are used as characteristic factors to perform multiple Logit regressions to establish a mapping relationship between traffic environment and categories. See Figure 4 for the implementation process. The implementation steps are as follows:

531) According to the method described in 51), a traffic characteristic model is established, and parameter estimation is performed.

 532) According to the previous finite mixed model parameter estimation result, the probability density function P, (x) of the travel speed distribution corresponding to the spatiotemporal sub-region on different dates is written, and the parameters are mixed Gaussian model as an example:

 κ

Pi ( ) =∑ j · ί] (ν _ξ ; Mj , _ί )

533) Calculate the Jensen-Shannon divergence d, _r between the two distributions.

 1

 d„ = JSD(P II Q) = -D(P \\ M) + -D(Q \\ M)

 twenty two

 Where Ρ, δ are two different probability distributions, Μ =— ΟΡ + β) , which is Kullback-Leibler divergence:

„ i ( )i.g In the case of a finite mixing model, the value cannot be explicitly expressed, but Monte Carlo sampling method can be used to approximate the calculation. The calculation method is:

534) Write the divergence between the two distributions as a distance matrix:

The matrix satisfies; = 4, / ⁼⁰ ( ^; '= ).

 535) Using the distance matrix as input to the K-Medoids algorithm, the clustering results are obtained, and the categories are indexed.

536) Using the category index as a response variable, the traffic environment data (including temperature, precipitation, visibility, etc.) is used as an independent variable. Perform multiple logit regression to obtain the mapping relationship between traffic environment E and traffic situation category T→7).

 537) Aggregate the same type of data, and re-establish the hybrid model with the new data set after aggregation, and perform parameter estimation to obtain the final historical traffic characteristic data set. The step 6) may specifically adopt the following methods:

 61) Simple real-time data processing. This method is implemented simultaneously with 51). The real-time traffic data is modeled and parameter estimated to obtain the characteristic function of the current traffic condition. The implementation steps of the method are exactly the same as 51), except that the data used is real-time traffic data.

 62) Classification processing. This method is carried out simultaneously with 52) or 53). Obtain the characteristic function of the traffic condition, and obtain the current information such as temperature, precipitation, visibility, traffic control measures, etc., and judge the current traffic condition category. See Figure 5 for the implementation process. The implementation steps are as follows:

621) calculating the travel speed in the time-space sub-zone, forming the overall _rt of the real-time travel speed _;

622) Establish a travel speed probability distribution model i^ (v _f , _rf ) = | - ί^, - μ^σ^, and perform parameter estimation;

623) Using the current traffic environment data (including temperature, precipitation, visibility, etc.) as an input parameter, using the mapping relationship ?0E→7 to obtain the category τ of the current traffic situation. The step 7) specifically includes the following steps:

 71) When step 62) is adopted, the historical traffic characteristic data under the category is located according to the category of the current traffic situation;

72) Calculate the difference between the two velocity distributions according to the description parameters η _Λ , μ _Λ , ^ of the current traffic characteristics and the description parameters η, μ, σ of the historical traffic characteristics: ^[(τ^ , μ^ σ^ ^η, μ, σ)^ ·/^^^ . When the historical traffic characteristics and real-time traffic characteristics (ie, the historical travel speed distribution and the real-time travel speed distribution) are similar, a smaller Jensen-Shannon divergence value will be obtained, that is, the difference between the two is smaller; when the historical traffic characteristics and When the real-time traffic characteristics are different, a larger Jensen-Shannon divergence value will be obtained, that is, the difference between the two is large, that is, the probability of existence of an abnormality is large, see FIG. The step 8) specifically includes the following steps:

 81) Normalize the difference in velocity distribution of each space-time sub-region to a normalized value of 0~1

 Diff^ - ra {diff)

 ξ' max diff, - min [diff,

82) Calculate the traffic anomaly index ^^ = >< 10 for each time and space sub-region. The present invention has the following advantages over similar technologies in the same field:

 (1) Make full use of the existing floating vehicle operation data (GNSS trajectory data), detect historical traffic state changes through historical traffic feature extraction and real-time traffic situation analysis, and realize real-time, low-cost, intelligent urban road traffic anomaly events. Detection

(2) Taking the probability distribution of the travel speed as the description of the traffic characteristics, the characteristics reflected are more comprehensive, avoiding the one-sidedness and instability of the traffic characteristics using a single index, and the reliability of detection is higher;

(3) Clustering a multi-logit regression combination for the characteristics of traffic characteristics affected by traffic environment (such as weather conditions) The algorithm establishes a mapping relationship between traffic environment characteristics and traffic situation categories;

 (4) According to the test of actual data, the urban road traffic anomaly detection technology based on floating car data proposed by the present invention can realize the detection of abnormal events with high accuracy, the detection rate exceeds 90%, and the false alarm rate is lower than 20%. It has achieved good detection results and can be applied to intelligent management and service of urban traffic. DRAWINGS

 The details and advantages of the present invention will become apparent and readily understood in conjunction with the following drawings in which:

 Figure 1 shows a schematic diagram of the components and basic principles of the present invention;

 Figure 2 is a schematic view showing the overall flow of the present invention in the implementation process;

 3 is a schematic diagram showing an implementation manner of a fast map matching algorithm of the present invention;

 4 is a schematic flow chart showing a historical traffic feature extraction scheme implemented by the present invention;

 FIG. 5 is a schematic flow chart showing a real-time traffic feature extraction scheme implemented by the present invention; FIG.

 Figure 6 shows a schematic diagram of the morphology of the Gaussian mixture model probability distribution;

 Figure 7 shows a measurement of the difference in the comparison between historical traffic characteristics and real-time traffic characteristics. Specific implementation

 In order to more clearly and clearly clarify the objects, technical solutions and advantages of the present invention, the specific embodiments of the present invention are described in detail below.

 As shown in Fig. 1, the overall system architecture of the present invention includes an onboard GNSS track recorder, a data center, a GNSS satellite, and a communication system carried by a floating vehicle. The GNSS here includes GPS, GLONASS, GALILEO, Beidou, IRNSS, QZSS and any similar navigation satellite positioning system. GNSS track recorders equipped with floating cars, buses, etc., with a certain sampling frequency / (general requirements; ).1Ηζ) record the position information of the vehicle at various points during driving, and through the GPRS mobile communication network (also The use of wireless network communication technologies such as WCDMA and TD-LTE, but the cost will be correspondingly improved), the location information is sent to the data center in real time. The data center establishes a historical road traffic characteristic database through data preprocessing, data fusion, and through a specific algorithm; establishes a real-time traffic feature database for the recently received real-time data; and determines whether the current traffic feature is abnormal through the mapping relationship between the historical database and the real-time database And visualize the display through the processing terminal and generate a traffic anomaly event report.

 The overall process of the scheme is shown in Figure 2, including the acquisition and storage of GNSS trajectory data, the establishment of spatiotemporal sub-areas, historical traffic feature extraction, real-time traffic feature extraction, and anomaly identification. Collecting and storing GNSS trajectory data is the data foundation of the whole scheme. Due to the huge amount of data, a distributed storage scheme should be adopted. For distributed storage, there are mature technologies, which are not the contents of the present invention. The basic assumption of establishing a spatiotemporal sub-area is that it has the same traffic characteristics in a specific area and a specific time period. This assumption is generally applicable after long-term observation. Historical traffic feature extraction, the principle is to use the GNSS trajectory data to calculate the travel speed, use a large number of travel speed data in the same space-time sub-region, establish a probability distribution model of vehicle speed, and estimate the parameters, and characterize the traffic characteristics with a small number of parameters. Real-time traffic feature extraction, the principle is to process and analyze the speed data in the current time period, and also establish the current vehicle speed probability distribution model. The anomaly identification is to use the difference measure to judge the degree of change of the real-time feature compared to the historical feature, and determine whether a traffic anomaly event occurs according to whether it reaches the threshold.

 According to the combination of the embodiments of the invention, the implementation is given below. Embodiment 1

Step 11. Determine the segment size of the time dimension by using the equidistant space-time division method, and the time segment span is a fixed value, usually 30 minutes. As a time segment; determine the segment size of the spatial dimension, the spatial segment span is a fixed value, usually takes a spatial grid of 200mX200m as a spatial segment.

 Step 12. Perform data preprocessing to perform data cleaning, data integration, data conversion, and data reduction on the GNSS positioning data to improve the structural degree of the data.

Step 13. Divide the spatial area to be processed into a grid of a certain size, and the range of each grid area can be expressed as = {(x _s , y _s ) I x _s e [x _r , x _r+l ) , y _s ≡ [y _r , y _r+l )}; determine the grid area where the anchor point is located, and use the distance and azimuth to search for the section where the anchor point is located; search for the section closest to point A, when the point A is satisfied If the difference between the driving direction angle and the direction angle of the road segment ij is less than the threshold value, that is, | - .|< is satisfied, the matching is completed; if the - -||< is not satisfied, the road segment zy is deleted in the search space, and the search for other road segments is continued. Until the condition is met; use the straight line equation of the road segment (if it is a curved road segment, it is roughly split into straight lines), calculate the projection coordinates of the GNSS positioning point on the road segment, and reduce the error caused by the GNSS positioning drift. The specific method is:

Determine the straight line equation of the road segment zy (if the road segment is a curve, divide it into several straight line segments): yy _i =k(xx _i ) where the slope is:

Xj - ^x , ) The projection line equation is: = _"( ^x _ ) ky _A — ky _j + k x. + x _A

 Solve the projected coordinates P as:

k ² y _A + y _t + h _A - h _t

' k ² +\

 After the map matching process, the anchor point is matched to the spatio-temporal sub-area in combination with the timestamp data of the coordinates of the positioning point.

Step 14. The overall vehicle speed data of each of the secondary floating cars in a time and space sub-region constitutes an overall. Calculate the speed of each car in the time-space sub-zone: = ^{2 + 3 +} ... + - , where ^... is the first and second of the time-space sub-zone

The distance between the GNSS positioning points, ..., the distance between the -1 and the nth GNSS positioning point, which is the first in the space-time sub-region, ..., the first Timestamp of the GNSS anchor point; the data in each spatio-temporal sub-area is not filtered to form a set for subsequent processing.

Step 15. As a whole, the historical data under the condition of no traffic anomaly is used to establish a traffic feature model and estimate the parameters. The method uses a finite mixing model to establish a traffic feature model and perform parameter estimation. Take the maximum component quantity K=3, and estimate the parameters of the mixed Gaussian model of «=1,2,... components separately; for f models, determine the best model by Bayesian information criterion (.BIO. Calculation: f ) =∑1⁄2f, ) n {\,2,...,5} A total of 5 mixed models. At the same time, calculate the five models BIC = -2\nL + k-\nn

 Where, is the maximum likelihood function value, t is the number of parameters in the model, "for the total amount of data.

After that, the /C smallest hybrid model is selected, and its parameter vectors 1], μ, and σ are recorded as the feature records of the present time-space sub-region. Step 16. Perform real-time traffic data model establishment and parameter estimation to obtain a characteristic function of the current traffic condition. The method is the same as step 1-5, and the parameter vectors t] _rt , μ<·ί, a _{rt are} recorded.

Step 17. Calculate the difference between the two velocity distributions according to the description parameters η _Λ , ^ of the current traffic characteristics and the description parameters 1], μ, σ of the historical traffic characteristics: ί3⁄4Γ[(η , μ _Λ , σ ), (η, μ, σ)] = JSDGP _rt ||P).

Step 18: normalize the speed of each space-time sub-region by ~1

 Calculate the traffic anomaly index ^^ = x 10 for each time-space sub-region. Embodiment 2

 Step 21: Using the equidistant space-time division method, determining the segment size of the time dimension, the time segment span is a fixed value, usually taking 30 minutes as a time segment; determining the segment size of the spatial dimension, the spatial segment span is a fixed value, usually taking 200m×200m The spatial grid acts as a spatial fragment.

 Step 22. Perform data preprocessing to perform data cleaning, data integration, data conversion, and data reduction on the GNSS positioning data to improve the structural degree of the data.

Step 23: Divide the spatial area to be processed into a grid of a certain size, and the range of each grid area can be expressed as = {(x _s , y _s ) I x _s e [x _r , x _r+l ) , y _s ≡ [y _r , y _r+l )}; determine the grid area where the anchor point is located, and use the distance and azimuth to search for the section where the anchor point is located; search for the section closest to point A, when the point A is satisfied If the difference between the driving direction angle and the direction angle of the road segment ij is less than the threshold value, that is, | - .|< is satisfied, the matching is completed; if the - -||< is not satisfied, the road segment zy is deleted in the search space, and the search for other road segments is continued. Until the condition is met; use the straight line equation of the road segment (if it is a curved road segment, it is roughly split into straight lines), calculate the projection coordinates of the GNSS positioning point on the road segment, and reduce the error caused by the GNSS positioning drift. The specific method is:

 Determine the straight line equation of the road segment zy (if the road segment is a curve, divide it into several straight line segments):

y. =k(xx.) where the slope is:

The projection line equation is: ^ =--τ( ^χ -3⁄4) ky _A ― ky _t + k ² x _t + x _A

 Solve the projected coordinates ρ

k ² y _A + y _t + h _A - h _t After the map matching process, the anchor point is matched to the spatio-temporal sub-area in combination with the timestamp data of the coordinates of the positioning point.

 d + d -d

 Step 24: Calculate the travel speed of each vehicle in the space-time sub-region: ― . where Α, 2... ί -1, " is the distance between the first and second GNSS positioning points in the space-time sub-area, . ....., the distance between the -1 and the GNSS anchor points, ^^ is the first time in the space-time sub-area, ..., the time stamp of the GNSS anchor points ; Specify the length of the time segment at the same time. The upper limit of the number of segments of the data segment ρ Search for the time data in a space-time sub-region. The velocity data in each time segment. If the number of velocity data in the time segment exceeds the upper limit ρ, the random data is added to V.

 Step 25. As a whole, the historical data under the condition of no traffic anomaly is used to establish the traffic characteristic model and estimate the parameters. The method uses a finite mixing model to establish a traffic feature model and perform parameter estimation. Take the maximum component quantity K=3 and estimate the parameters of the mixed Gaussian model of «=1,2,... components separately; for f models, determine the best model by Bayesian information criterion (.BIO. Calculation: f ) =∑1⁄2f, ) n {\,2,...,5} A total of 5 mixed models. At the same time, calculate the five models

 Where, is the maximum likelihood function value, t is the number of parameters in the model, "for the total amount of data.

 After that, the /C smallest hybrid model is selected, and its parameter vectors 1], μ, and σ are recorded as the feature records of the local space-time sub-region. According to the parameter estimation result, the probability density function p'Q) of the travel speed distribution corresponding to the spatio-temporal sub-region on different dates is written:

Calculate the Jensen-Shannon divergence between the two pairs of distributions:

 1

 4 = JSD(P II Q) = -D(P \\M) + -D(Q \\M)

 twenty two

 Where Ρ, δ are two different probability distributions, M = OP + (3⁄4, for Kullback-Leibler divergence:

D{P\\Q) = ±P{x _k )\o _g ^{ In the case of a finite-mixing model, Monte Carlo sampling is used to approximate the calculation. The calculation method is:

 Write the divergence between the two distributions into a distance matrix:

 d d.

 D- d d This matrix satisfies ί3⁄4=φ, d, corpse Q = j,.

The distance matrix is used as the input of the K-Medoids algorithm to obtain clustering results and index the categories. Using the category index as the response variable, the traffic environment data (including temperature, precipitation, visibility, etc.) is used as an independent variable, and multiple logit regressions are performed to obtain the mapping relationship R(E→T) between the traffic environment E and the traffic situation category T.

 The same type of data is aggregated, and the hybrid model is re-established with the new data set after aggregation, and the parameter estimation is performed to obtain the final historical traffic characteristic data set.

 Step 26: Obtain a characteristic function of the traffic condition, and obtain current information such as temperature, precipitation, visibility, traffic control measures, and the type of the current traffic condition.

Calculate the travel speed in the time-space sub-zone, constitute the real-time travel speed overall V _irt; establish the travel speed probability distribution model pJ _{Vi rt} ) = ^η . · f people ν _ξ , _μι , σ , and make parameter estimation; The data (including temperature, precipitation, visibility, etc.) is used as an input parameter to obtain the category T of the current traffic situation using the mapping relationship R(E→T).

Step 27: Locating the historical traffic characteristic data of the category according to the current traffic situation category ;; calculating the two speed distributions according to the description parameters η, μ, σ of the current traffic characteristics description parameters 3⁄4, ^, ^ and historical traffic characteristics The difference: diff [(η _Λ , μ _Λ , ) , (η, μ, σ)] = JSD(P _rt \\ Ρ). Step 28: Denormalize the speed of each time and space sub-region

 Calculate the traffic anomaly index ^^ = x 10 for each time-space sub-region. Embodiment 3

Step 31: Using a non-equidistant space-time division method, for a central area of the city where the road network density is greater than ² km/km ² or the peak hour flow is greater than 1000 vehicles/hour, take a 30 min time segment and a 200 mX 200 m spatial segment for the road network density. For suburban areas of less than 2km/km ² or peak hour traffic of less than 1000 vehicles per hour, 30 min time segments and 400 mX 400 m space segments are taken.

 Step 32: Perform data preprocessing, perform data cleaning, data integration, data conversion, and data reduction on the GNSS positioning data to improve the structural degree of the data.

Step 33: Divide the spatial area to be processed into a grid of a certain size, and the range of each grid area may be expressed as = {< , ) I [ ^x _r ^ ^x _r+1 ly _s [n)};

The floating vehicle GNSS data acquisition frequency is expressed as f ₀ =\lt ₀ , and the point adjacent to A in time (ί4-.). :) is defined as 1-adjacent point of A, ■P ko), /3⁄44+23⁄4) is defined as 2 - neighboring point of A, and so on, then corpse ((4-), (4+3⁄4)) is defined as The neighboring point of A. in/. When <1Ηζ, take Λ=1 or 2. Take the road segment ij with the smallest distance from the neighboring point of the distance Α and Α, and calculate the mean value of the driving direction angle of the 1-adjacent point of A and A ^ 4, if |^4. - | < is satisfied, complete the match; otherwise, search for other road segments Until the conditions are met.

 Use the straight line equation of the road segment (if it is a curved road segment, it is roughly split into straight lines), calculate the projection coordinates of the GNSS positioning point on the road segment, and reduce the error caused by the GNSS positioning drift. The specific method is:

Determine the straight line equation of the road segment zy (if the road segment is a curve, divide it into several straight line segments): y - y = k{x - x The slope is:

投影直线方程为: y-y_A τ(^χ_ ) ky_A - ky_t + k²x_t + x_A

The projection line equation is: yy _A τ( ^χ _ ) ky _A - ky _t + k ² x _t + x _A

 Solve the projected coordinates Ρ:

k ² y _A + y _t + h _A - h _t

y _P

 + i

 After the map matching process, the anchor point is matched to the spatio-temporal sub-area in combination with the timestamp data of the coordinates of the positioning point.

ά, _Ί +ά _Ί

 Step 34: Calculate the travel speed of each vehicle in the space-time sub-region: , where <3⁄4, 2...ί -1," is the distance between the first and second GNSS positioning points in the space-time sub-region, . ....., the distance between the -1 and the GNSS anchor points, ... is the first time in the space-time sub-region, ..., the time of the GNSS anchor point Poke; specify the time segment length at the same time segment number of data segments upper limit ρ search for a time-space sub-region time ζ·speed data within each time segment, if the number of velocity data within the time segment exceeds the upper limit ρ random data is added to V .

 Step 35: Perform historical traffic data without traffic abnormality as a whole, and establish traffic feature model and parameter estimation. The method uses a finite mixing model to establish a traffic feature model and perform parameter estimation. Take the largest number of components K = 3 and estimate the parameters of the mixed Gaussian model of «=1, 2, ... components separately; for f models, determine the best model by Bayesian information criterion (.BIO). :

/ ^ν ) = ί /=1 Λ( ^ν ) "ε{ΐ,2"..,5} A total of 5 mixed models. At the same time, calculate the five models

 BIC = -2\nL + k-\nn

 Where, is the maximum likelihood function value, t is the number of parameters in the model, "for the total amount of data.

 After that, the /C smallest hybrid model is selected, and its parameter vectors 1], μ, and σ are recorded as the feature records of the local space-time sub-region. According to the parameter estimation result, the probability density function p'Q) of the travel speed distribution corresponding to the spatio-temporal sub-region on different dates is written:

A ( ) = ∑ Λ ( ) Calculate the Jensen-Shannon divergence d, _r between the two distributions.

d _u = JSD(P II Q) = -D(P \\M) + -D(Q \\M)

 twenty two

Where Ρ, δ are two different probability distributions, Μ = - OP + ρ) , which is Kullback-Leibler divergence:

In the case of a finite mixing model, Monte Carlo sampling is used to approximate the calculation. The calculation method is: 3⁄4 _c ( ll g) = -∑ig Writes the divergence between the two pairs as a distance matrix:

 The matrix satisfies ί3⁄4=φ, d, corpus Q = j,.

 The distance matrix is used as the input of the K-Medoids algorithm to obtain clustering results and index the categories.

 Using the category index as the response variable, the traffic environment data (including temperature, precipitation, visibility, etc.) is used as an independent variable to perform multiple logit regression to obtain the mapping relationship between the traffic environment E and the traffic situation category T R{E→T).

 The same type of data is aggregated, and the hybrid model is re-established with the new data set after aggregation, and the parameter estimation is performed to obtain the final historical traffic characteristic data set.

 Step 36: Obtain a characteristic function of the traffic condition, and obtain current information such as temperature, precipitation, visibility, traffic control measures, and the type of the current traffic condition.

Calculate the travel speed in the time-space sub-zone, which constitutes the real-time travel speed overall V _irt; establish the travel speed probability distribution model _Prt ( _{V rt} ) = ^η . · f people ν _ξ , _μι , σ , and make parameter estimation; The environmental data (including temperature, precipitation, visibility, etc.) is used as an input parameter to obtain the category T of the current traffic situation using the mapping relationship R(E→T).

Step 37: Locating the historical traffic characteristic data of the category according to the current traffic situation category ;; calculating the two speed distributions according to the description parameters of the current traffic characteristics 3⁄4, ^, ^ and the description parameters η, μ, σ of the historical traffic characteristics The difference: diff [(η _Λ , μ _Λ , ) , (η, μ, σ)] = JSD(P _rt \\ Ρ). Step 38: Denormalize the speed of each time and space sub-region

 Calculate the traffic anomaly index ^^ = x 10 for each time-space sub-region.

Claims

Claim  1. An urban road traffic anomaly detection method based on floating car data, comprising the following steps:  1) Establish a spatio-temporal sub-area: divide the day into a number of time segments, and divide the implementation area of urban road traffic anomaly detection into several spatial segments;  2) Preprocessing of historical trajectory data: processing the floating vehicle GNSS positioning history data into the sampled vehicle speed data of the historical trajectory;  Preprocessing of real-time trajectory data: processing floating vehicle GNSS positioning real-time data into sampled vehicle speed data of real-time trajectory;  3) Historical trajectory data analysis and feature extraction: Using the sampled vehicle speed data of the historical trajectory, establishing a historical travel speed probability distribution, and obtaining a historical traffic characteristic model;  Real-time trajectory data analysis and feature extraction: using the sampled vehicle speed data of the real-time trajectory, establishing a real-time travel speed probability distribution, and obtaining a real-time traffic feature model;  4) Anomaly detection: The difference between the historical traffic characteristic model and the real-time traffic characteristics is measured by Jensen-Shannon divergence, and the difference between historical and real-time traffic characteristics is obtained;  5) Quantitative characterization of abnormal severity: Using the historical and real-time traffic feature differences, the traffic condition anomaly index is calculated. 2. The method for detecting urban road traffic anomaly based on floating car data according to claim 1, wherein the step 1) Use one of the following methods:  La) Isometric space-time division method: Determine the segment scale of the time dimension, the time segment span is a fixed value, usually take 30mm as a time segment; determine the segment scale of the spatial dimension, the spatial segment span is a fixed value, usually takes 200m X 200m The spatial grid acts as a spatial segment; Lb) Non-equidistant space-time division method: For urban central areas with road network density greater than 2km/km ² or peak hour flow greater than 1000 vehicles/hour, take 30mm time segments and 200m X200m space segments for road network density less than 2km /km ² or a suburb of a city with a peak hour flow of less than 1000 vehicles per hour, taking a 30mm time segment and a 400m X 400m space segment. 3. The method for detecting urban road traffic anomaly based on floating car data according to claim 1, wherein the step 2) The preprocessing of the historical trajectory data includes:  2a) Data structuring: data cleaning, data integration, data conversion, and data reduction of the GNSS positioning history data of the floating vehicle to obtain structured GNSS positioning history data;  2b) Fast map matching: Combine the urban road network data, map the structured GNSS positioning historical data to the urban road network through the map matching algorithm, and establish the matching relationship between the positioning points and the road segments in the structured GNSS positioning historical data. Obtaining a matching relationship table between the positioning point and the road segment in the structured GNSS positioning history data, and correcting an error caused by the positioning drift;  2c) Vehicle speed calculation and sampling: Calculate the vehicle speed based on the structured GNSS positioning history data, obtain historical vehicle speed data, and perform data sampling on the historical vehicle speed data to obtain sampled historical vehicle speed data. Replacement page (Article 26) The method for detecting an urban road traffic anomaly based on floating car data according to claim 1, wherein the preprocessing of the real-time trajectory data in step 2) comprises:  2d) Data structuring: Data processing, data integration, data conversion, and data reduction of GNSS positioning real-time data of floating vehicles are obtained, and structured GNSS positioning real-time data is obtained;  2e) Fast map matching: Combine the urban road network data, map the structured GNSS positioning real-time data to the urban road network through the map matching algorithm, and establish the matching relationship between the positioning points and the road segments in the structured GNSS positioning real-time data. Obtaining a matching relationship table between the positioning point and the road segment in the structured GNSS positioning real-time data, and correcting an error caused by the positioning drift;  2f) Vehicle speed calculation and sampling: Calculate the vehicle speed based on the structured GNSS positioning real-time data, obtain real-time vehicle speed data, and perform data sampling on the real-time vehicle speed data to obtain sampled real-time vehicle speed data. The method for detecting an urban road traffic anomaly based on floating car data according to claim 3, wherein the fast map matching in step 2b) comprises:  2b 1) Divide the space area to be processed into a grid of a certain size, and the range of each grid area can be expressed as 2b2) Determine the grid area where the anchor point is located, and use the distance and azimuth to search for the road segment where the anchor point is located; 2b3) Calculate the projection coordinates of the GNSS anchor point on the road segment by using the GNSS anchor point linear projection method. The method for detecting an urban road traffic anomaly based on floating car data according to claim 4, wherein the fast map matching in step 2e) comprises:  2e l) Divide the space area to be processed into a grid of a certain size, and the range of each grid area can be expressed as 4 = {( ^x ^ I [x ₊ i), y, [υ ₊₁ ]}  2e2) Determine the grid area where the anchor point is located, and use the distance and azimuth to search for the road segment where the anchor point is located; 2e3) Calculate the projection coordinates of the GNSS anchor point on the road segment by using the GNSS anchor point linear projection method. 7. The method for detecting urban road traffic anomaly based on floating car data according to claim 5, wherein step 2b2) adopts one of the following methods:  2b21) Single point matching method: Search for the nearest road segment from point A. When the difference between the traveling direction angle satisfying point A and the direction angle of the road segment ij is less than the threshold value, then | - | < is satisfied, and the matching is completed; if not satisfied | - | < , delete the road segment in the search space and continue to search for other road segments until the conditions are met. 2b22) Point sequence matching method: This scheme is applicable to high frequency floating car data. The floating vehicle GNSS data acquisition frequency is expressed as /ϋ=1>, and the point PO^o), Pfe+ o) adjacent to the time is defined as the 1-adjacent point of the ,, P(t _A -2h), Ρθ4 +2. ) is defined as the 2-adjacent point of A, and so on, then Pfe-feo), Pi + kto) is defined as the / 1-adjacent point of A; <When 1Ηζ, fetch / fe = l or 2; / 1 adjacent the point of minimum distance from the segment A and A ij is taken, and calculates the average of A ^ and A in the direction of travel adjacent the corner points', if yes _- θ υ \ < δ _θ , complete the match; otherwise, search for other segments until the condition is met. Replacement page (Article 26) The method for detecting urban road traffic anomaly based on floating car data according to claim 3, wherein the preprocessing of the historical trajectory data in step 2c) adopts one of the following methods:  2cl) Full sample method: The overall speed data of each sub-floating vehicle in a time-space sub-region constitutes the total. The implementation method is to calculate the travel speed of each vehicle in the space-time sub-zone: = w + . + ''", Where 2. . . ^^ is the distance between the 1st and 2nd GNSS anchor points in the space-time sub-region, ..., the distance between the -1 and the GNSS anchor points , ... „ is the first in the time and space sub-area, ..., the first Timestamp of the GNSS anchor point; the data in each spatio-temporal sub-area is not filtered to form a set for subsequent processing; 2c2) Time-smooth sampling method: Specify the length of the time segment and the upper limit of the number of segments of the same time; Search for the velocity data in each time segment in a time-space sub-region. If the number of velocity data in the time segment exceeds the upper limit, the upper limit is randomly selected. The data of the number of bars is added to the data sample to be processed. The implementation method is to calculate the travel speed of each vehicle in the space-time sub-zone: = W -""- ¹ '", where ₂ ...^, „ is the time zone The distance between the 1st and the 2nd GNSS positioning point, ..., the distance between the n-1th and the nth GNSS positioning point, ... is the first in the space-time sub-area, .. ..., the "timestamp of the GNSS anchor point; the specified time segment length at the same time segment data strip number upper limit search for speed data in the time segment of a time-space sub-region, if the speed data within the time segment Exceeding the upper limit; 1⁄2 _∞ , random data is added! ^. The method for detecting an urban road traffic anomaly based on floating car data according to claim 4, wherein the preprocessing of the historical trajectory data in step 2f) adopts one of the following methods: 2fl) Full sample method: The total vehicle speed data of each sub-floating vehicle in a time-space sub-region constitutes the whole. The implementation method is to calculate the travel speed of each vehicle in the space-time sub-zone: = W ..""- ¹ '" , among them 4, ₂ . The distance between the 1st and 2nd GNSS positioning points in the space-time sub-region, ..., the distance between the -1 and the GNSS positioning points, ... „The first in the time and space sub-area, ..., the first Timestamp of the GNSS anchor point; the data in each spatio-temporal sub-area is not filtered to form a set for subsequent processing; 2f2) Time-smooth sampling method: Specify the length of the time segment and the upper limit of the number of segments of the same time; Search for the velocity data in each time segment in a time-space sub-region. If the number of velocity data in the time segment exceeds the upper limit, the upper limit is randomly selected. The data of the number of bars is added to the data sample to be processed. The implementation method is to calculate the travel speed of each vehicle in the space-time sub-region: V^ W -""- ¹ '", where 2...^, „ is the time-space sub-region The distance between the first and second GNSS anchor points, the distance between the -1 and the nth GNSS anchor points, Replacement page (Article 26) ... is the first time in the time-space sub-area, ..., the time stamp of the GNSS anchor point; the upper limit of the number of segments of the specified time segment length at the same time; 1⁄2 _∞ ; searching for a space-time sub-region The speed data in the time segment of the time, if the number of speed data in the time segment exceeds the upper limit; 1⁄2 _∞ , the random data is added to ^. 10. The method for detecting urban road traffic anomaly based on floating car data according to claim 1, wherein the historical trajectory data analysis and feature extraction in step 3) adopt one of the following methods: 3a) Simple historical trajectory data fusion method: The historical data under the condition of no traffic anomaly, as a whole, the traffic feature model establishment and parameter estimation, the method uses the finite hybrid model to establish the traffic feature model and perform parameter estimation;  3b) Historical trajectory data classification method: According to temperature, precipitation, visibility and traffic control measures, the historical data under the condition of no traffic anomalies is divided into different categories, and models and parameter estimation are respectively established; 3c) Historical data clustering method: For historical data, through the comparison between time and space sub-regions, obtain the quantitative quantitative representation of different time and space regions, and use the quantized differences for clustering; temperature, precipitation, visibility and As a characteristic factor, traffic control measures carry out multiple Lo _gl regressions to establish a mapping relationship between traffic environment and categories. 11. The method according to claim 10, wherein the step 3a) adopts one of the following methods:  3al) Fixed-component mixed Gaussian model method: The mixed Gaussian model with a fixed number of components is used to describe the probability distribution of the vehicle speed. The number of components is manually specified according to the distribution pattern of the vehicle speed under typical conditions. In order to ensure the reliability of the probability distribution, the number of components cannot be too small. , generally available = 4~6;  3a2) Mixed Gaussian model method with variable number of components: The method based on model evaluation is used to select the appropriate number of components as follows: Determine the maximum number of possible components, and separately mix the components of “=1, 2,...^ The Gaussian model performs parameter estimation; for each model, the best model is determined by Bayesian information criteria; 3a3) Finite mixed model method in which the number of components and the type of distribution are variable: the distribution pattern of the subcomponents and the number of components are variable. The method for detecting an urban road traffic anomaly based on floating car data according to claim 10, wherein the historical data clustering method according to step 3c) comprises:  3d) performing parameter estimation of the historical traffic characteristic model;  3c2) according to the parameter estimation result, write a probability density function p, (x) of the travel speed distribution corresponding to the spatio-temporal sub-region on different dates; 3c3) Calculate the Jensen-Shannon divergence φ _μ between the two distributions  = SD{P II Q) = ^D(P II M) + ^ D(Q \\ M) where \ β is two different probability distributions, Μ = 1(Ρ + 0, for Kullback-Leibler Degree: Replacement page (Article 26) „ i ( )i.g In the case of a finite mixing model, the value cannot be explicitly expressed, but Monte Carlo sampling method can be used to approximate the calculation. The calculation method is: D _MC (f II g) = -∑log ^{( i)} ― D(f \\ g)  n g 3c4) Write the divergence between the two distributions as a distance matrix:

 The matrix satisfies 4= , =0 (/=/); 3c5) Using the distance matrix as input to the K-Medoids algorithm, obtaining clustering results, and indexing the categories; 3c6) Using the category index as the response variable, using traffic environment data (including temperature, precipitation, visibility, etc.) as independent variables , perform multiple Lo _gl t regression, obtain the mapping relationship between traffic environment E and traffic situation category T (E→T) _; 3c7) aggregate data of the same category, and re-establish the hybrid model by using the new data set after aggregation And perform parameter estimation to obtain the final historical traffic characteristic data set. The method for detecting urban road traffic anomaly based on floating car data according to claim 10, wherein the real-time trajectory data analysis and feature extraction in step 3) adopt one of the following methods:  3d) Simple real-time data processing method: This method is implemented simultaneously with the simple historical trajectory data fusion method described in 3a), and the real-time traffic data is modeled and parameter estimated to obtain the characteristic function of the current traffic condition;  3e) Classification processing: This method is implemented simultaneously with the historical trajectory data classification method described in 3b) or the historical data clustering method described in 3c), obtaining the characteristic function of the traffic condition, and simultaneously obtaining the current temperature and precipitation. Information such as quantity, visibility, traffic control measures, and the type of current traffic conditions. The method for detecting an urban road traffic anomaly based on floating car data according to claim 13, wherein the classification processing method according to step 3e) comprises: 3el) Calculate the travel speed in the time-space sub-zone, which constitutes the overall _rt of the real-time travel speed _; 3e2) Establish a travel speed probability distribution model p _rt (ν _{ξ Γ} , ) = | · f ( _r ,; μ , ), and perform parameter estimation; 3e3) Using the current traffic environment data (including temperature, precipitation, visibility, etc.) as input parameters, use the mapping relationship ?(£→7) to obtain the category of the current traffic situation. The method for detecting an urban road traffic anomaly based on floating car data according to claim 1, wherein: step 4) the abnormality detecting comprises:  4a) According to the category T of the current traffic situation, locate the historical traffic characteristic data under the category, if there is no category, Replacement page (Article 26) It is not necessary to distinguish categories; 4b) Calculate the difference between the two velocity distributions according to the description parameters τ^, μ„, o _{rt of the} current traffic characteristics and the description parameters η, μ, σ of the historical traffic characteristics: ί3⁄4Γ[(η , μ , σ ), (η ,μ, _σ )] = \\Ρ). 16. The method for detecting urban road traffic anomaly based on floating car data according to claim 1, wherein: step 5) quantitative characterization of abnormal severity includes:  5a) Normalize the difference in velocity distribution for each spatiotemporal sub-region to a normalized value of 0~1:  Diff^ -ram{diff)  Max {diff, - min {diff,  5b) Calculate the traffic anomaly index for each time-space subzone = > 10. Replacement page (Article 26)