A method for obtaining high temporal and spatial resolution port traffic observation data from AIS data

Abstract

The present invention discloses a method for obtaining high-temporal and spatial resolution port traffic observation data from AIS data, comprising the following steps: S1, preprocessing historical ship trajectories based on AIS data; S2, dividing the preprocessed historical ship trajectories into dwell segments and movement segments using a dwell index formula and a K-means algorithm; extracting the center point coordinates of the dwell segments using an isolation forest algorithm to obtain ship travel chains; labeling each dwell segment with its dwell type and the port to which it belongs based on spatial rules and process logic; S3, constructing key indicators reflecting port traffic status through temporal and spatial statistics of berthing and anchoring behaviors; constructing a directed shipping network with ports as nodes, and using network analysis methods to explore traffic connections and characteristics between ports. This method can improve the timeliness and accuracy of port traffic situation monitoring, overcoming the problems of slow update, coarse granularity, and high heterogeneity of traditional statistical data.

Description

Technical Field

[0001] The present invention relates to the technical field of port traffic, and in particular to a method for obtaining high-temporal and spatial resolution port traffic observation data from AIS data. Background Art

[0002] Since 2004, the International Maritime Organization (IMO) has mandated that all ships over 300 gross tonnage be equipped with AIS (Automatic Identification System) equipment. This allows shipping companies and port authorities to obtain multi-dimensional real-time data, including location, speed, and status. Unlike traditional methods that rely on annual port reports or static statistics, AIS data is highly frequent and continuous, making it suitable for analyzing granular issues such as vessel dynamics, traffic distribution, and path evolution. It has been widely used in research areas such as port service efficiency assessment, ship emissions, and maritime network resilience. Using AIS data, port traffic data can be constructed from a ship's perspective, enabling optimized routes, increased turnover rates, and improved channel design and berth utilization.

[0003] Mining and analyzing ship activity patterns from AIS data using a data-driven approach has become a research hotspot in the field of transportation engineering. Existing research mainly includes three categories: grid-based, statistical, and vector-based. Grid-based methods divide the traffic area into discrete grids and analyze ship traffic characteristics within each grid. While suitable for identifying traffic hotspots, this method is susceptible to grid division when analyzing large areas and is computationally expensive. Statistical methods use statistical methods to analyze the distribution of ship traffic attributes. While computationally efficient, they typically only analyze the distribution of a few variables, making it difficult to uncover complex behaviors hidden within trajectories and difficult to visualize. Vector-based methods treat trajectories as path vectors, analyzing traffic flows and typical routes. Chinese patent application publication number CN117785963A discloses a method and system for mining maritime traffic networks based on open-source AIS data. This method can accurately extract traffic networks over large sea areas from AIS data. It offers advantages such as high extraction accuracy, low manual supervision costs, and rich semantic patterns, enhancing maritime situational awareness and understanding. However, this method places high demands on the extraction of feature points, such as origin and destination points (OD points), and many studies fail to consider the specific traffic meaning and contextual information of these points, hindering a deeper understanding of ship behavior and traffic patterns.

[0004] In summary, although the above three methods have their own advantages, they also have limitations. The existing methods fail to fully utilize AIS data to accurately extract ship activities between ports. In trajectory modeling, there is also a general lack of in-depth description of behavior, making it difficult to combine attribute information in the trajectory. Summary of the Invention

[0005] Purpose of the invention: The purpose of the present invention is to provide a method for obtaining high-temporal and spatial resolution port traffic observation data from AIS data. The method uses AIS data to construct terminal-level ship OD data, records the minute-level arrival time and berthing time of ships, extracts the ship travel chain, realizes the automatic recognition of berthing, anchoring and other behaviors and port ownership judgment, and improves the timeliness and accuracy of port operation status identification.

[0006] Technical solution: A method for obtaining high-temporal and spatial resolution port traffic observation data from AIS data, comprising the following steps:

[0007] S1, based on AIS data, performs denoising, interpolation and time resampling operations on the ship's historical trajectory in sequence to obtain the preprocessed ship's historical trajectory;

[0008] S2: Based on the stay index formula and K-means algorithm, the pre-processed ship historical trajectory is divided into stay segments and movement segments. The center point coordinates of the stay segments are extracted using the isolation forest algorithm to obtain the ship travel chain. Based on spatial rules and process logic, each stay segment is labeled with the stay type and the port to which it belongs.

[0009] S3, through the time-space statistics of berthing and anchoring behaviors, constructs key indicators reflecting the port traffic status to describe the port congestion situation and operation efficiency; constructs a directed shipping network with ports as nodes, and combines network analysis methods to explore the traffic connections and characteristics between ports.

[0010] Furthermore, in step S2, the stay index formula is used to calculate the stay probability of the trajectory point, with a value range of 0 to 1. A larger value indicates a higher likelihood of stay. The expression is as follows:

[0011] ,

[0012] in, For trajectory points The stay index; Is the size of the sliding window, indicating the trajectory points; For trajectory points and The distance between is the attenuation coefficient of spatial distance; exp() represents the natural exponential function;

[0013] Next, the K-means algorithm was applied to cluster the stay indices of all trajectory points. The initial cluster center was set to [0, 1], and the trajectory was divided into "movement segments" and "stay segments". Adjacent stay segments were merged based on time continuity to obtain complete stay segments. The rule for merging stay segments is as follows: for a movement segment between two adjacent stay segments, if the duration of the movement segment is less than 60 minutes and the maximum distance is less than 1 km, it is judged as an invalid movement segment and merged.

[0014] Furthermore, in step S2, the isolation forest algorithm is used to eliminate abnormal points in the stop segment and extract the center point coordinates of the stop segment; then, the start and end time, center point coordinates, and corresponding MMSI and ship attribute fields of each stop segment are extracted to form a structured "ship travel chain" dataset.

[0015] Furthermore, based on the port information and ship attributes, each stopover segment is marked with the stopover type and the port to which it belongs. The specific implementation steps are as follows:

[0016] SA1: Draw the berth boundary based on the remote sensing image and the clustering features of the center points of the stop segment. Points falling within the berth boundary are identified as "berth stop segments". The specific implementation steps are as follows:

[0017] SA11: Identify the clustered areas of the center points of the stopover segments, combine remote sensing images to determine whether they are close to the coastline and whether they have typical port facilities, draw the terminal boundaries, and form a vector surface;

[0018] SA12, using ArcGIS tools to construct buffer zones for the dock boundaries, setting the inland river buffer zone to 50 meters and the coastal buffer zone to 100 meters;

[0019] SA13, extract the center points of the stay segments that fall into the buffer zone through the spatial screening method and mark them as berth stays;

[0020] SA2: Based on the identified anchorage stops, the non-anchorage segments preceding each anchorage stop segment are extracted as candidate anchorage stop segments. The isolation forest algorithm is used to perform an anomaly score on each candidate anchorage stop segment. The lower the score, the more likely it is non-anchorage behavior. The true anchorage stop segments are then screened based on the threshold.

[0021] SA3 uses the city-level administrative divisions to determine the port to which each berth stay belongs. For an anchorage stay, if the next operation section of a certain anchorage stay is a berth stay within a port, the anchorage stay is determined to belong to the same port.

[0022] Furthermore, in step S3, key indicators reflecting the port traffic status are constructed, including:

[0023] The geographic hash area is used to measure the spatial range covered by the center points of all stopover segments within the statistical time unit. By statistically analyzing the number of all unique hash codes within the time period, the effective spatial coverage area of ​​the port during the statistical time period is estimated. It is defined as follows:

[0024] ,

[0025] in, Indicates the number of unique GeoHash codes, Represents the area of ​​each grid;

[0026] Gross ton-hour load is used to assess the port's capacity load per unit time. For each stopover segment, the proportion of its stopover time within the statistical time period is calculated as a weighting coefficient, which is then multiplied by the deadweight tonnage of the ship. The total load of all stopover ships is then summed to obtain the overall load of the port area during the statistical time period. The calculation formula is as follows:

[0027] ,

[0028] in, For the The length of time a ship stays during the statistical period, is the length of the statistical time period, For the The deadweight tonnage of the vessel; Indicates the total number of ships in the statistical time period;

[0029] The number of ships arriving, staying, and leaving the port is used to reflect the dynamic traffic status of the port within a specific time period. Based on the arrival time, departure time, and stay type of each stop segment, the number of arrivals, departures, and stays are counted separately. The calculation formula is as follows:

[0030] ,

[0031] ,

[0032] ,

[0033] in, is the statistical time period; Indicates the number of arrivals, Indicates the number of departures, Indicates the number of ports; Count() is a counting function used to count the number of records that meet the conditions;

[0034] The average turnaround time is used to assess the efficiency of ships' operations in port. The turnaround time of each ship is calculated from the start and end times of its stay, and then divided into the corresponding statistical time period according to the arrival time. The final result is the average value of all ships in the statistical time period. The calculation formula is as follows:

[0035] ,

[0036] in, Indicates the The actual length of time the ship stayed in the port.

[0037] Furthermore, in step S3, the directed shipping network is used to describe the interconnection structure between ports. Each edge in the network represents a connection where a ship actually sails. The weight of the edge is defined as the cumulative number of ships that actually sail from port p to port q, which represents the freight intensity of the directed route. By counting the cumulative number of voyages between port pairs, a port OD voyage matrix is ​​formed, in which the elements represents the cumulative number of voyages from port p to port q; network analysis methods are used to quantify network indicators from three perspectives: connection strength, network centrality index and community detection.

[0038] Compared with the prior art, the present invention has the following significant effects:

[0039] 1. This paper constructs terminal-level ship OD data based on raw AIS data, recording the minute-by-minute arrival time and berthing duration of ships, and extracting the ship travel chain. First, it combines preprocessing steps such as denoising, interpolation, and resampling to significantly improve data quality. Then, it integrates the stay index, cluster analysis, and anomaly detection to achieve automatic recognition of berthing and anchoring behaviors and port ownership judgment, enhancing the structured expression and semantic understanding capabilities of trajectory data.

[0040] 2. This paper constructs a port operation index system covering dimensions such as operating pressure, traffic density, and spatial distribution. Based on the index data, it establishes a port transportation network and uses network analysis methods to reveal the connection structure and operating characteristics between ports at multiple levels, enriching the analytical dimensions of port transportation data.

[0041] 3. The present invention obtains port traffic observation data with high temporal and spatial resolution by constructing four key indicators that can reflect the port traffic status and conducting statistical analysis. This can improve the timeliness and accuracy of port operation status identification, overcome the problems of slow updating, coarse granularity and strong heterogeneity of traditional statistical data, and provide reliable support for smart port construction and shipping big data applications. BRIEF DESCRIPTION OF THE DRAWINGS

[0042] Figure 1 is a flow chart of the present invention;

[0043] Figure 2 It is the flow chart of the dwell segment identification algorithm;

[0044] Figure 3 This is a schematic diagram of the calculation of the stay index based on the sliding window;

[0045] Figure 4 It is a line graph of the total number of stay segments and the average total stay duration under different parameter combinations, where (a) shows the change in the total number of stay segments under different parameter combinations, and (b) shows the change in the average total stay duration under different parameter combinations;

[0046] Figure 5 This is a schematic diagram of the effect of the isolation forest algorithm, where (a) shows the points eliminated by the algorithm in a single berthing and anchoring trajectory, and (b) shows the comparison of the space of stopover points in Nanjing Port and Shanghai Port before and after the spatial streamlining in December 2023;

[0047] Figure 6 It is a flowchart of the berth identification algorithm and an example diagram of its application to two types of terminals (coastal / inland river);

[0048] Figure 7 This is the outlier identification effect of Suzhou Port anchor points based on GMM isolation forest score and the KDE graph of location and time interval variables;

[0049] Figure 8 This is the outlier identification effect of Nanjing Port anchorage points based on GMM isolation forest score and the KDE diagram of location and time interval variables;

[0050] Figure 9 This is the outlier identification effect of Shanghai Port anchorage points based on GMM isolation forest score and the KDE graph of location and time interval variables;

[0051] Figure 10 This is the outlier identification effect of Yingkou Port anchorage points based on GMM isolation forest score and the KDE diagram of location and time interval variables;

[0052] Figure 11 This is a multi-dimensional statistical chart of the traffic operation characteristics of national container ports in 2023, where (a) is a trend chart of the average daily operation indicators of national ports, and (b) is a comparison chart of the average monthly turnover time of ports in various regions;

[0053] Figure 12 This is a normalized comparison chart of national container port centrality indicators in 2023;

[0054] Figure 13 This is the spatial distribution map of the connectivity strength among container ports across China in 2023;

[0055] Figure 14 This is the result of identifying the community structure of national container ports in 2023;

[0056] Figure 15 This is a relationship diagram of container transport flows among various container port communities across the country in 2023. DETAILED DESCRIPTION

[0057] The present invention will be described in further detail below with reference to the accompanying drawings and specific implementations.

[0058] like Figure 1 A method for obtaining high temporal and spatial resolution port traffic observation data from AIS data comprises the following steps:

[0059] Step 1: Reconstruct the historical trajectory of the ship. First, perform denoising, interpolation, and time resampling on the AIS data to obtain the preprocessed historical trajectory of the ship.

[0060] Step 11, trajectory denoising;

[0061] This step aims to remove drift points in the trajectory data, that is, parts of the trajectory that show spikes or abnormal protrusions, thereby improving the overall quality of the ship's trajectory data. AIS data records the instantaneous state information of the ship at each moment. The speed and heading relative to the ground carried by it often have errors and are difficult to use directly to identify abnormal trajectory points. Therefore, the present invention first uses the spatial position and timestamp of the trajectory points to calculate the speed (SOG) and rate of turn (ROT) between two adjacent points, and generates a set of data describing the dynamic characteristics of each ship for each ship. Subsequently, the box plot method is used to determine the abnormal detection thresholds of the ship's speed and turning rate, respectively, as follows:

[0062] (1)

[0063] (2)

[0064] Among them, T is the maximum threshold of each variable, Q1 and Q3 are the 25% quantile and 75% quantile of a set of data, and IQR is the interquartile range.

[0065] By setting the SOG and ROT thresholds, and applying a "segmentation-filtering" strategy to remove outliers, the system removes outliers. Specifically, when both the speed and turning rate of a trajectory segment exceed their respective thresholds, it is identified as an outlier. Breakpoints are then set between the start and end points of the segment, dividing the entire trajectory into multiple sub-trajectories. The system then checks the number of track points in each sub-trajectory, removing those with fewer than 10 points and concatenating the remaining sub-trajectories to complete the trajectory denoising process.

[0066] Step 12, trajectory interpolation;

[0067] By judging the distance between adjacent trajectory points, we can determine whether interpolation is needed to improve the continuity and spatiotemporal resolution of the trajectory. The specific method is as follows: when the distance S between two adjacent trajectory points is greater than 2 kilometers, it is determined that there is a large gap in the section and interpolation is required. First, determine the timestamp set of the required interpolation points , which is calculated as follows:

[0068] (3)

[0069] (4)

[0070] (5)

[0071] Among them, the collection Indicates the timestamp of the interpolation point, ; 、 are the timestamps of the starting and ending points of the missing segment, S is the straight-line distance between the starting and ending points of the missing segment, and St is the preset distance interval threshold, which is 1 km.

[0072] After determining the timestamp, the coordinates of each interpolation point are calculated using the Lagrange linear interpolation method. Let the starting point of the missing segment be , the end point is , the timestamp corresponding to the interpolation point is , then the coordinates corresponding to the interpolation points are It can be given by formula (6) and (7):

[0073] (6)

[0074] (7)

[0075] Step 13, time resampling;

[0076] Furthermore, because the uneven temporal sampling frequency of the original AIS data hinders subsequent trajectory analysis and behavior recognition, the present invention resamples the interpolated trajectory data at one-minute intervals. This process not only helps compress the data volume but, more importantly, aligns the time intervals between trajectory points, facilitating the subsequent calculation of features such as velocity and acceleration, and simplifying related formulas and processing flows.

[0077] After systematic data cleaning and preprocessing in step 1, high-quality historical ship trajectory data with strong spatiotemporal continuity was obtained. However, this data still appears as a sequence of points, lacking clear behavioral labels and semantic orientation, exhibiting the typical characteristics of "high spatiotemporal and low semantics." The various behavioral characteristics exhibited by ships during navigation, such as departure, underway, berthing, anchoring, and waiting, cannot be directly captured from historical trajectories. Therefore, it is urgent to introduce trajectory mining and enhance semantic information to achieve structured analysis of trajectory data.

[0078] Step 2: Divide the pre-processed ship's historical trajectory into stop segments and movement segments; extract the arrival time, departure time, and center location of each stop segment to construct a ship travel chain; and combine port geographic information to annotate each stop segment with the stop type (berth stop / anchorage stop) and the port to which it belongs. Step 21: Extract the stop segments from the trajectory using the stop index formula and the K-means algorithm;

[0079] This step proposes a “stay index” to measure the stay degree of a trajectory in its life cycle based on the spatial clustering characteristics of trajectory points, and constructs an automatic stay segment identification method based on this. The flowchart is as follows: Figure 2 shown.

[0080] Specifically, the sliding window method is first used to calculate the stay index of each trajectory point, with the window range being the n trajectory points before and after. Since step 13 has resampled the trajectory to a fixed 1-minute interval, the calculations of different trajectory points are comparable. The stay index formula is used to calculate the stay probability of the trajectory point. The stay index range is [0,1]. The higher the value, the greater the probability that the point is a stay. The schematic diagram of calculating the stay index of each trajectory point is shown in the figure below. Figure 3 As shown, the calculation formula is as follows:

[0081] (8)

[0082] in, For trajectory points The dwell index; 2n+1 is the size of the sliding window, representing the n trajectory points before and after; For trajectory points and The distance between is the attenuation coefficient of spatial distance; exp() represents the natural exponential function.

[0083] Formula (8) follows the first law of geography of spatial influence attenuation and embodies the idea of ​​kernel density analysis. Figure 4 ), the optimal parameter values ​​are determined to be the sliding window size n=20 and the decay coefficient =60, this set of parameters shows the best effect.

[0084] Next, the K-means algorithm was applied to cluster the dwell indices of all trajectory points, with the initial cluster centers set to [0, 1]. Trajectories were divided into "movement segments" and "dwelling segments." To avoid misidentification due to short movement segments, adjacent dwell segments were merged based on temporal continuity to construct complete dwell segments. The specific merging rules are as follows: For a movement segment between two adjacent dwell segments, if the duration of the movement segment is less than 60 minutes and the maximum distance is less than 1 kilometer, it is considered an invalid movement segment and is merged.

[0085] The automatic identification method of the stop segment can effectively identify the ship's stop behaviors such as anchoring and berthing. Figure 2 As shown in the figure, taking the trajectory of a ship in December 2023 as an example, the K-means algorithm identified 13 stops and 7 voyages. The results show that the automatic stop segment identification method has good accuracy and practicality in identifying stops and dividing voyage segments.

[0086] Step 22: Using the isolation forest algorithm, the coordinates of the center points of the stop segments are extracted to obtain the ship travel chain.

[0087] The stop segments extracted in step 21 often contain some transition points and noise points, such as trailing, spurring, or drifting points at the end of the trajectory. If these outliers are not processed and directly used in the center point calculation, they will cause spatial position offsets and reduce the spatial accuracy of the stop segments. To address this issue, the present invention introduces an isolation forest algorithm to remove spatial outliers within each stop segment, retaining only representative core trajectory points. The center coordinates are then calculated based on these points, which more accurately reflect the ship's actual stop location. This invention implements the isolation forest algorithm using the IsolationForest class provided by the scikit-learn Python library and sets the outlier ratio parameter (contamination) to 0.2, identifying and removing the least representative 20% of trajectory points in each stop segment. Finally, the start and end times, center point coordinates, and corresponding MMSI (Maritime Mobile Service Identity) and vessel attribute fields of each stop segment are extracted to form a structured "ship travel chain" dataset for subsequent port traffic data statistics.

[0088] Figure 5 The comparison of the effects before and after the application of the isolation forest algorithm is shown. Figure 5 There are obvious outliers in the original berthing and anchoring sections shown in (a). After processing, the center point positioning is more accurate. Figure 5 Panel (b) further shows the distribution of center points of all stopover sections at Nanjing Port and Shanghai Port in December 2023. After screening using the isolation forest algorithm, the center points are concentrated in typical stopover areas, significantly improving the spatial differentiation between anchorage and berthing.

[0089] Step 23, adding semantic information, marking the type of stay and the port to which it belongs for each stay segment;

[0090] To enhance the semantic recognition of ship stop chains, the present invention, building on the stop segment extraction in step 21, further labels each stop segment with its stop type (berth stop / anchorage stop) and associated port. The key to achieving this goal lies in constructing berth and anchorage boundaries, then classifying stops based on their center locations: those falling within the berth boundary are designated as "berthing segments," those falling within the anchorage boundary are designated as "anchoring segments," and the remainder are treated as noise and discarded. Berth boundaries are drawn based on the clustering features of remote sensing imagery and stop segment centers, providing a clear geographic location. Anchorage boundaries are inferred based on a spatiotemporal model, representing the reasonable range and time window from the anchorage to the berth, assisting in identifying the spatial characteristics of anchoring behavior. Port identification relies on the terminal's Area of ​​Interest (AOI) attribution. By correlating the processes between berthing and anchoring, each stop segment is accurately assigned to a specific port.

[0091] Step 231, berth stay extraction;

[0092] First, the terminal AOI areas of all research ports are extracted, and the adjacent stay segments are screened based on the terminal AOI areas. See the flowchart for details. Figure 6 Specifically, it includes: taking the center point of each stop section as the representative location, visually overlaying it on the map, identifying the spatial clustering area of ​​the center points of the stop section, and combining remote sensing images to determine whether it is close to the coastline and whether it has typical port facilities (such as yards, quay cranes, etc.) to confirm the berth location; then drawing the berth boundary based on the clustered area to form a vector surface; then using ArcGIS tools to construct a buffer zone for the berth boundary, considering different water conditions, setting the inland river buffer zone to 50 meters and the coastal zone to 100 meters; finally, using spatial screening methods to extract the center points of the stop section that fall into the buffer zone, and marking the stop section as a berth stop.

[0093] Step 232, anchorage stay extraction;

[0094] The core of anchorage stay identification is to restore the typical port entry process of "anchorage waiting-berth operation". Under normal circumstances, the effective anchorage stay segment appears before the berth stay segment, and assumes the function of queuing and waiting for operation. Therefore, based on the identified berth stay segment, the adjacent non-berth stay segment in each trajectory is traced back as a candidate anchorage stay segment. However, not all stop behaviors before berth stay belong to anchoring. Some may be short stops, or the ship may enter the port directly without actual waiting. Therefore, the present invention introduces the isolation forest algorithm to identify outliers, and combines the constraint theory in time geography to construct four characteristic variables including spatial position and time distance:

[0095] Table 1 Input features of anchor recognition algorithm

[0096]

[0097] First, for each port, candidate anchorage stay segments corresponding to berths of all terminals in the port are screened, and the above four features (i.e., the longitude of the center point of the stay segment, the latitude of the center point of the stay segment, the time interval from the end of the current stay segment to the beginning of the next segment, and the straight-line distance from the center point of the current stay segment to the center point of the terminal) are calculated. Next, the isolation forest algorithm is used to perform an abnormality score on each candidate anchorage stay segment. The lower the score, the more likely it is an atypical anchoring behavior. However, since the abnormal proportion is unknown and a fixed threshold cannot be preset, the present invention introduces a Gaussian mixture model (GMM). Based on the property that the isolation forest score approximates a normal distribution, multiple Gaussian distributions are adaptively fitted, and the mean plus the standard deviation of the last distribution ( ) as the discrimination threshold. The GMM model requires a preset number of distributions. Considering that most terminals usually correspond to one anchorage, their score distribution should be unimodal. However, in port areas such as Shanghai Port, where anchorages are dense and ship types are diverse, there may be multiple anchorages and corresponding multiple distribution peaks. To adapt to such complexity, the present invention automatically selects the optimal distribution number from candidate numbers 1 to 6 through the Bayesian Information Criterion (BIC), and determines the final distribution accordingly. As the outlier threshold, all candidate anchorage segments with scores lower than the threshold are finally deleted.

[0098] Figure 7 、 Figure 8 、 Figure 9 、 Figure 10 The results of the present invention on identifying abnormalities at candidate anchorages in Suzhou Port, Nanjing Port, Shanghai Port, and Yingkou Port are shown. The first row shows the isolation forest score distribution and GMM fitting results. GMM adaptively determines the optimal number of peaks and adds the mean plus standard deviation of the last peak ( ) is used to set a threshold. The second row shows the KDE distribution plot (Kernel Density Estimation) of the location feature and the anomaly identification results. Anchorage points often correspond to high-density areas, and some ports have multiple anchorage areas. The number of peaks is consistent with the GMM fitting results. The third row shows the KDE plot and identification results of the time interval feature. Multiple time clustering areas also correspond to multiple Gaussian peaks, and points with significant time deviations are effectively marked as outliers.

[0099] Step 233, the port to which it belongs is determined;

[0100] This method combines administrative boundary information with city-level administrative divisions to determine the port to which each terminal AOI belongs. For berth stays, since their center points have been filtered through a spatial buffer and clearly fall within the terminal AOI of a specific port, they can be directly assigned to that port. For anchorage stays, this method analyzes the actual ship operation chain and subsequent associated berth stays to determine anchorage segment ownership. Specifically, if the next operation segment of a certain anchorage stay is a berth stay within a port, the anchorage stay is determined to belong to the same port. This method not only ensures accurate port identification but also reflects the consistency of operational logic within the port.

[0101] Step 3: statistics of port traffic observation data;

[0102] Step 31, port operation status;

[0103] This paper constructs an indicator system for port traffic monitoring, aiming to comprehensively characterize port operations from both spatial and temporal dimensions. This indicator system is based on the spatial and temporal statistics of ship berthing and anchoring behavior, and proposes four core indicators to describe port congestion and operational efficiency.

[0104] 1) Geohash Area: This measures the spatial extent covered by the center points of all stopover segments within a statistical time unit. Geohash, a spatial encoding method proposed by Gustavo Niemeyer, features good hierarchical and adjacency properties. In this example, a coding level with a precision of 7 is used, with each grid covering approximately 0.34 km² of sea area. By statistically analyzing the number of unique hash codes within a time period, the effective spatial coverage area of ​​the port during the statistical period can be estimated, defined as follows:

[0105] (9)

[0106] in, Indicates the number of unique GeoHash codes, =0.34 km² for each grid area.

[0107] 2) Gross Ton-Hour Load: This is used to assess the port's capacity load per unit time. Specifically, for each stopover segment, the proportion of its duration within the statistical period is calculated as a weighting factor, which is then multiplied by the ship's deadweight tonnage (DWT). The total is then calculated for all stopover ships to obtain the overall port load during the statistical period:

[0108] (10)

[0109] in, For the The length of time a ship stays during the statistical period, is the length of the statistical time period, is the deadweight tonnage of the ship, and M represents the total number of ships in the statistical time period.

[0110] 3) Number of vessels arriving, staying, and leaving: This reflects the dynamic traffic status of the port within a specific timeframe. Based on the arrival time, departure time, and stay type of each stop segment, statistics are compiled for three dimensions: arrivals (by start time), departures (by end time), and stays (both start and end times cover the statistical timeframe). By categorizing the stay types, traffic dynamics can also be calculated for berths and anchorages:

[0111] (11)

[0112] (12)

[0113] (13)

[0114] in, is the statistical time period; Indicates the number of arrivals, Indicates the number of departures, Indicates the number of ports; Count() is a counting function used to count the number of records that meet the conditions.

[0115] 4) Average Turnaround Time: This is used to assess the efficiency of a ship's operations in port. Each ship's turnaround time is calculated from the start and end times of its stay, and then allocated to the corresponding statistical time period based on its arrival time. The final result is the average value of all ships in the statistical time period:

[0116] (14)

[0117] in, Indicates the The actual length of time the ship stayed in the port. Indicates the The arrival time of this ship, Indicates the The departure time of the ship, M represents the total number of ships in the statistical time period.

[0118] The average turnaround time can be calculated at the berth, anchorage and overall port area levels to enable efficiency comparison and port operation monitoring.

[0119] The above four indicators can be statistically analyzed and visualized at multiple temporal and spatial scales. Figure 11 As shown, Figure 11 (a) shows the changing trends of the number of daily arriving ships, average turnaround time, gross tonnage-hour load, and GeoHash area across China; Figure 11 (b) shows the radar chart of the monthly average berth, anchorage and total port area turnover time of inland ports and coastal ports.

[0120] Step 32, port shipping network situation;

[0121] Furthermore, the present invention constructs a directed shipping network with ports as nodes, which is used to systematically describe the connectivity structure and freight flow relationship between ports. Each node in the network represents a port, and the directed edges between nodes indicate the connection relationship of actual ship navigation. The weight of the edge is defined as the cumulative number of ships actually sailing from port p to port q, which represents the freight intensity of the directed route. To construct the network, the "berth stay segment" information of all ships is first extracted from the stay chain, and the two consecutive berth stays of each ship are paired in chronological order to form the actual transportation segment between the starting and ending ports. Each segment record represents a real navigation event of a ship from one port to another. By summarizing the segment information of all ships and counting the cumulative number of voyages between each port pair, a port pair (Origin-Destination, OD) voyage matrix is ​​generated. Each element of the matrix represents the cumulative number of voyages from port p to port q.

[0122] Subsequently, network analysis methods were applied to quantify the structural relationships between ports along three dimensions: first, the strength of connectivity was calculated, including the port's in-degree (indicating the number of connected ports) and weighted in-degree (indicating the total number of port voyages); second, network centrality metrics were assessed, such as closeness centrality (measuring port accessibility), betweenness centrality (demonstrating its importance in transit routes), and eigenvector centrality (reflecting its overall influence); third, the community detection algorithm (Louvain) was used to identify cluster structures within the network, revealing potential regional logistics circles and functional zoning between ports. This analysis method was implemented using open-source network analysis tools such as NetworkX and Community in Python. The community detection module utilized the python-louvain package, while centrality metrics and connectivity strength calculations were based on relevant functions provided by NetworkX.

[0123] Figures 12 to 15 Shows the results of network analysis from multiple angles: Figure 12 The normalized scores of each major port under four centrality indicators (including degree centrality, weighted degree, betweenness centrality and eigenvector centrality) can be used to identify the core nodes in the network; Figure 13 The spatial distribution pattern of the connection intensity between ports is shown, revealing high-intensity transportation corridors; Figure 14 To use the Louvain algorithm to detect communities on the network, multiple geographically closely related regional clusters are formed, which has strong spatial explanatory power; Figure 15 The use of chord diagrams to display the intensity of ship transportation between communities helps to identify the flow interaction relationship between the regional collaborative transportation system and the core ports.

[0124] With the above-mentioned ideal embodiment of the present invention as inspiration, through the above description, relevant personnel can make various changes and modifications without departing from the scope of the technical idea of ​​this invention.

Claims

1. A method for obtaining high temporal and spatial resolution port traffic observation data from AIS data, characterized in that: The following steps are involved: S1, based on AIS data, performs denoising, interpolation and time resampling operations on the ship's historical trajectory in sequence to obtain the preprocessed ship's historical trajectory; S2: Based on the stay index formula and K-means algorithm, the pre-processed ship historical trajectory is divided into stay segments and movement segments. The center point coordinates of the stay segments are extracted using the isolation forest algorithm to obtain the ship travel chain. Based on spatial rules and process logic, each stay segment is labeled with the stay type and the port to which it belongs. S3: By analyzing the temporal and spatial statistics of berthing and anchoring behaviors, key indicators reflecting port traffic status are constructed to describe port congestion and operational efficiency. A directed shipping network with ports as nodes is constructed, and network analysis methods are used to explore the traffic connections and characteristics between ports. In step S2, the stay index formula is used to calculate the stay probability of the trajectory point, with a value range of 0 to 1. A larger value indicates a higher likelihood of stay. The expression is as follows: ， in, For trajectory points The stay index; Is the size of the sliding window, indicating the trajectory points; For trajectory points and The distance between is the attenuation coefficient of spatial distance; exp() represents the natural exponential function; Next, the K-means algorithm was applied to cluster the stay indices of all trajectory points, with the initial cluster center set to [0,1], dividing the trajectory into "movement segments" and "stay segments." Adjacent stay segments were merged based on temporal continuity to obtain complete stay segments. The rule for merging stay segments is as follows: for a movement segment between two adjacent stay segments, if the movement segment lasts less than 60 minutes and the longest distance is less than 1 kilometer, it is considered an invalid movement segment and merged. In step S3, key indicators reflecting the port traffic status are constructed, including: The geographic hash area is used to measure the spatial range covered by the center points of all stopover segments within the statistical time unit. By statistically analyzing the number of all unique hash codes within the time period, the effective spatial coverage area of ​​the port during the statistical time period is estimated. It is defined as follows: ， in, Indicates the number of unique GeoHash codes, Represents the area of ​​each grid; Gross ton-hour load is used to assess the port's capacity load per unit time. For each stopover segment, the proportion of its stopover time within the statistical time period is calculated as a weighting coefficient, which is then multiplied by the deadweight tonnage of the ship. The total load of all stopover ships is then summed to obtain the overall load of the port area during the statistical time period. The calculation formula is as follows: ， in, For the The length of time a ship stays during the statistical period, is the length of the statistical time period, For the The deadweight tonnage of the vessel; Indicates the total number of ships in the statistical time period; The number of ships arriving, staying, and leaving the port is used to reflect the dynamic traffic status of the port within a specific time period. Based on the arrival time, departure time, and stay type of each stop segment, the number of arrivals, departures, and stays are counted separately. The calculation formula is as follows: ， ， ， in, is the statistical time period; Indicates the number of arrivals, Indicates the number of departures, Indicates the number of ports; Count() is a counting function used to count the number of records that meet the conditions; The average turnaround time is used to assess the efficiency of ships' operations in port. The turnaround time of each ship is calculated from the start and end times of its stay, and then divided into the corresponding statistical time period according to the arrival time. The final result is the average value of all ships in the statistical time period. The calculation formula is as follows: ， in, Indicates the The actual length of time the ship stayed in the port.

2. The method for obtaining high temporal and spatial resolution port traffic observation data from AIS data according to claim 1, characterized in that: In step S2, the isolation forest algorithm is used to remove outliers within the stop segments and extract the center point coordinates of the stop segments. Then, the start and end times, center point coordinates, and corresponding MMSI and ship attribute fields of each stop segment are extracted to form a structured "ship travel chain" dataset.

3. The method for obtaining high temporal and spatial resolution port traffic observation data from AIS data according to claim 1, characterized in that: According to the port information and ship attributes, each stopover segment is marked with the stopover type and the port to which it belongs. The specific implementation steps are as follows: SA1: Draw the berth boundary based on the remote sensing image and the clustering features of the center points of the stop segment. Points falling within the berth boundary are identified as "berth stop segments". The specific implementation steps are as follows: SA11: Identify the clustered areas of the center points of the stopover segments, combine remote sensing images to determine whether they are close to the coastline and whether they have typical port facilities, draw the terminal boundaries, and form a vector surface; SA12, using ArcGIS tools to construct buffer zones for the dock boundaries, setting the inland river buffer zone to 50 meters and the coastal buffer zone to 100 meters; SA13, extract the center points of the stay segments that fall into the buffer zone through the spatial screening method and mark them as berth stays; SA2: Based on the identified anchorage stops, the non-anchorage segments preceding each anchorage stop segment are extracted as candidate anchorage stop segments. The isolation forest algorithm is used to perform an anomaly score on each candidate anchorage stop segment. The lower the score, the more likely it is non-anchorage behavior. The true anchorage stop segments are then screened based on the threshold. SA3 uses the city-level administrative divisions to determine the port to which each berth stay belongs. For an anchorage stay, if the next operation section of a certain anchorage stay is a berth stay within a port, the anchorage stay is determined to belong to the same port.

4. The method for obtaining high temporal and spatial resolution port traffic observation data from AIS data according to claim 1, characterized in that: In step S3, the directed shipping network is used to describe the connectivity structure between ports. Each edge in the network represents a connection where a ship actually sails. The weight of the edge is defined as the cumulative number of ships that actually sail from port p to port q, which represents the freight intensity of the directed route. By counting the cumulative number of voyages between port pairs, a port OD voyage matrix is ​​formed, in which the elements represents the cumulative number of voyages from port p to port q; network analysis methods are used to quantify network indicators from three perspectives: connection strength, network centrality index, and community detection.

Images