A human-machine collaboration and confrontation control optimization method for complex dynamic environment

Abstract

The application discloses a kind of man-machine cooperation and confrontation control optimization methods for complex dynamic environment, it is related to man-machine collaborative control technical field.The method obtains machine observation data, human behavior data and machine control output data, constructs corresponding time sequence and feature vector sequence;According to sampling time, sampling period and transmission delay, determine causal reachable relationship and construct causal reachable relationship matrix, establish man-machine closed-loop benchmark path model based on historical benchmark data;Predict the current state using benchmark path parameters, obtain residual vector sequence and determine confrontation event anchor point;According to the difference of event before and after cross-source residual correlation matrix, extract event after mixed features, and obtain causal correction state time series data through path attribution.The application can reduce the influence of asynchronous data mismatch and event after mixed on control decision, improve the accuracy and stability of man-machine cooperation and confrontation control in complex dynamic environment.

Description

Technical Field

[0001] This invention relates to the field of human-machine collaboration technology, and more specifically to a human-machine collaboration and adversarial control optimization method for complex dynamic environments. Background Technology

[0002] With the development of intelligent robots, engineering equipment, and unmanned systems, human-machine collaborative control is gradually evolving from single manual operation or autonomous control to human-environment collaborative decision-making. In complex dynamic environments, systems typically receive machine-observed data such as environmental conditions and equipment motion states collected by machine-side sensors, as well as data generated by operators through handle inputs, gestures, electromyography, or other behavioral data. This data is then combined with machine control outputs to continuously adjust the actuators. Because machine-side data usually has a high sampling frequency, while human behavior is influenced by perception, cognition, and action response processes, its sampling period and response delay differ significantly from machine-side data, making it difficult to establish a stable one-to-one time correspondence between data from different sources.

[0003] Existing human-machine collaborative control methods typically employ timestamp matching, resampling, or interpolation to synchronize data from different sources, and then use the synchronized data for state recognition, intent judgment, or control strategy calculation. However, in complex dynamic environments, whether data have the same or similar sampling times does not directly reflect their actual sequential relationship. Communication delays, device response times, and human cognitive reaction times can all cause dynamic temporal shifts between different data. Directly associating data with adjacent times can easily lead to using data without actual interaction as the basis for control, affecting the accurate judgment of human-machine state relationships.

[0004] Furthermore, in the closed-loop operation of human-machine interaction, human actions can influence machine control and machine state, and the machine control results, in turn, further influence subsequent human behavior through environmental changes and operational feedback. When sudden obstacles, load changes, external disturbances, or sudden changes in the strategy of adversaries occur in the environment, human behavior, machine observation state, and machine control output may change simultaneously, generating new correlations that do not exist in the normal operation phase. Existing methods typically identify state changes based on pre-established fixed relationship models or single abnormal indicators, making it difficult to distinguish whether the changes are caused by human active operation, machine self-control feedback, or external adversarial disturbances. This can easily lead to misinterpreting common changes after an event as the original human-machine interaction relationship. When these abnormal correlations are directly incorporated into subsequent control strategy calculations, it may cause the control system to incorrectly attribute the human's true operational intent, the machine's autonomous control action, and external disturbances, resulting in inappropriate adjustments to control parameters. In highly dynamic human-machine collaboration scenarios, this problem may further manifest as frequent switching of control authority, repeated correction of operating instructions, oscillating machine responses, and control strategy drift, thereby reducing the stability and adaptability of human-machine collaboration and adversarial control in complex dynamic environments. Summary of the Invention

[0005] The purpose of this invention is to provide a method and system for recognizing the usage status of pet supplies based on image classification, thereby solving the problems existing in the prior art.

[0006] To achieve the above objectives, the present invention provides the following technical solution: a human-computer collaboration and adversarial control optimization method for complex dynamic environments, comprising: S1. Acquire machine observation data, human behavior data, and machine control output data; construct the first, second, and third time series according to the sampling time; and extract the first, second, and third feature vector sequences. S2. Determine the causal reachability relationship between feature vectors based on the sampling time, sampling period and transmission delay of each time series, and construct the causal reachability relationship matrix; S3. Based on the first, second and third eigenvectors that satisfy the causal reachability relationship within the historical benchmark time window, establish a human-machine closed-loop benchmark path model and obtain the benchmark path parameter matrix. S4. Predict the first and second eigenvectors within the current time window using the baseline path parameter matrix to obtain the first and second residual vector sequences, and determine the adversarial event anchor point based on their joint mutation. S5. Based on the adversarial event anchor points, construct the pre-event cross-source residual correlation matrix and the post-event cross-source residual correlation matrix respectively. Perform difference operation on the two to obtain the cross-source residual increment matrix, and extract the post-event hybrid features that reach the judgment threshold from it to form a hybrid feature vector sequence. S6. Based on the baseline path parameter matrix and the sequence of confounding feature vectors, determine the first, second and third state contribution components corresponding to human behavior, machine control output and post-event confounding factors, and determine the path attribution coefficient vector with state reconstruction error as a constraint. S7. Correct the contribution components of the first, second and third states based on the path attribution coefficient vector, generate causal correction state time series data, and input it into the control optimization module to generate machine control instructions to execute human-machine cooperation and adversarial control.

[0007] As can be seen from the above technical solution, the present invention has the following beneficial effects: This invention asynchronously processes machine observation data, human behavior data, and machine control output data using time-series feature processing. By establishing causal reachability relationships based on sampling time, sampling period, and transmission delay, it reduces erroneous associations caused by hard alignment of data at different time scales. Furthermore, it predicts the current human-machine closed-loop state using historical baseline paths and identifies new confounding features generated by adversarial events through changes in cross-source residuals before and after events. This distinguishes the impact of human behavior, machine control output, and external disturbances on machine state changes. Correcting the state contribution components based on path attribution results reduces the probability of external disturbances and closed-loop feedback being misjudged as genuine human intent or normal machine responses, generating causal-corrected state time-series data that reflects the current actual operating state. Using this causal-corrected state time-series data for control optimization improves the accuracy and stability of control decisions in complex dynamic environments, reduces the probability of human-machine control conflicts, frequent control switching, and closed-loop oscillations, thereby enhancing the adaptability of human-machine collaborative and adversarial control systems to sudden disturbances and dynamic environmental changes. Attached Figure Description

[0008] Figure 1 A flowchart of a human-computer collaboration and adversarial control optimization method for complex dynamic environments provided by the present invention; Figure 2 A timing diagram for a human-computer collaboration and adversarial control optimization method for complex dynamic environments provided by this invention; Figure 3 This is a graph showing the trajectory and joint mutation index of the present invention. Figure 4 This is a comparison chart of the RMSE of position tracking under different resistance intensities according to the present invention; Figure 5 This is a comparison chart of system recovery time after an adversarial event under different update strategies of the present invention; Figure 6 This is a schematic diagram of the cross-source residual increment matrix of the present invention. Detailed Implementation

[0009] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0010] Example 1: like Figure 1 and Figure 2 As shown, the present invention provides a technical solution: a human-computer collaboration and adversarial control optimization method for complex dynamic environments, comprising: S1. Acquire machine observation data, human behavior data, and machine control output data; construct the first, second, and third time series according to the sampling time; and extract the first, second, and third feature vector sequences. S2. Determine the causal reachability relationship between feature vectors based on the sampling time, sampling period and transmission delay of each time series, and construct the causal reachability relationship matrix; S3. Based on the first, second and third eigenvectors that satisfy the causal reachability relationship within the historical benchmark time window, establish a human-machine closed-loop benchmark path model and obtain the benchmark path parameter matrix. S4. Predict the first and second eigenvectors within the current time window using the baseline path parameter matrix to obtain the first and second residual vector sequences, and determine the adversarial event anchor point based on their joint mutation. S5. Based on the adversarial event anchor points, construct the pre-event cross-source residual correlation matrix and the post-event cross-source residual correlation matrix respectively. Perform difference operation on the two to obtain the cross-source residual increment matrix, and extract the post-event hybrid features that reach the judgment threshold from it to form a hybrid feature vector sequence. S6. Based on the baseline path parameter matrix and the sequence of confounding feature vectors, determine the first, second and third state contribution components corresponding to human behavior, machine control output and post-event confounding factors, and determine the path attribution coefficient vector with state reconstruction error as a constraint. S7. Correct the contribution components of the first, second and third states according to the path attribution coefficient vector, generate causal correction state time series data, and input it into the control optimization module to generate machine control instructions to execute human-machine cooperation and adversarial control. S8. From the self-adversarial event anchor point, stop using post-event data to update the human-machine closed-loop baseline path model and enter the parameter freeze state; The path credibility of the current data is calculated based on the path attribution coefficient vector and the state reconstruction error. When the time window for a consecutive preset number of joint mutation indicators is lower than the recovery threshold and the attribution coefficient corresponding to post-event confounding factors is lower than the preset limit, the parameter freeze state is lifted. Select samples whose path credibility reaches the update threshold from the data after the parameter freeze is lifted, and use them as credible update samples. The baseline path parameter matrix is ​​updated incrementally using trusted update samples to suppress the contamination of parameters in the human-machine closed-loop baseline path model by post-event mixed data.

[0011] In the above scheme, human-machine multi-source data acquisition, causal reachability determination, human-machine closed-loop baseline path modeling, adversarial event identification, post-event hybrid feature extraction, state contribution attribution, causal correction control, and baseline path model freezing and reliable update are completed sequentially according to S1 to S8. This separates the abnormal correlation data after the adversarial event from the normal human-machine closed-loop data, and generates machine control instructions using the causally corrected state data during the control process.

[0012] Specifically, in S1, three data sources are used: machine observation data, human behavior data, and machine control output data during machine operation. Machine observation data includes the machine's position, speed, acceleration, attitude, target distance, target relative speed, and sensor detection values. Human behavior data includes joystick displacement, direction input, button status, pedal input, and touch input. Machine control output data includes drive force commands, braking force commands, steering commands, actuator displacement commands, torque commands, and speed commands. The sampling period for machine observation data is set to 20 milliseconds, the sampling period for human behavior data is set to 50 milliseconds, and the sampling period for machine control output data is set to 20 milliseconds. Each sampled data point simultaneously records the data source number, sampling time, reception time, and sampled value.

[0013] Machine observation data are arranged in chronological order of sampling time to form the first time series; human behavior data are arranged in chronological order of sampling time to form the second time series; and machine control output data are arranged in chronological order of sampling time to form the third time series. For each data channel, the current sampled value, the rate of change of adjacent sampled values, and the arithmetic mean of five consecutive sampled values ​​are extracted. The rate of change of adjacent sampled values ​​is calculated by subtracting the previous sampled value from the current sampled value, and then dividing this change by the time difference between the current and previous sampling times. The arithmetic mean of five consecutive sampled values ​​is calculated by adding the current sampled value to the previous four consecutive sampled values ​​and then dividing by 5. The current sampled value, rate of change, and arithmetic mean are connected according to a pre-fixed data channel order to obtain the first, second, and third feature vectors, respectively. The feature vectors are arranged according to their corresponding sampling times to form the first, second, and third feature vector sequences, respectively.

[0014] In the above scheme, in S2, the original sampling times of the first, second, and third time series are retained, and the three time series are not forcibly converted to the same sampling interval through interpolation. The actual transmission delay corresponding to each sampled data is obtained by subtracting the data sampling time from the data reception time. For any source feature vector and any target feature vector, the sampling time of the source feature vector is first added to its actual transmission delay to obtain the earliest time when the source feature vector can reach the subsequent processing node. Then, this earliest arrival time is compared with the sampling time of the target feature vector. When the earliest arrival time of the source feature vector is not later than the sampling time of the target feature vector, and the time difference between the sampling time of the target feature vector and the sampling time of the source feature vector does not exceed 500 milliseconds, the source feature vector and the target feature vector are determined to be causally reachable. When either of the above two conditions is not met, the source feature vector and the target feature vector are determined to be causally unreachable.

[0015] Further, the corresponding discrete lag positions are determined according to the sampling period of the source time series. Specifically, the time difference between the sampling time of the target feature vector and the sampling time of the source feature vector is divided by the sampling period of the source time series; if the division result has a decimal, it is rounded up to the nearest integer, and the integer is used as the lag order of the source feature vector relative to the target feature vector. This judgment is then applied to all possible connections between the first, second, and third feature vectors. Matrix positions that satisfy the causal reachability condition are marked as 1, and matrix positions that do not satisfy the causal reachability condition are marked as 0. All judgment results form a causal reachability relation matrix. In subsequent human-machine closed-loop baseline path modeling, path parameters are only established for positions with a matrix value of 1, while the path parameters for positions with a matrix value of 0 are fixed at 0.

[0016] In the above scheme, during initial modeling in S3, the historical baseline time window is selected from the continuous 60 seconds of running data after the system starts up. After the system has been running and completed at least one model update, the historical baseline time window is composed of the continuous 60 seconds of data before the current time that has not been marked as adversarial event data.

[0017] For each target feature vector within the historical benchmark time window, the first, second, and third feature vectors within a 500-millisecond range are extracted, and features that do not satisfy the causal reachability relation with the target feature vector are deleted according to the causal reachability relation matrix. The retained features are concatenated according to the data source order, time lag order, and feature dimension order to form the historical input data corresponding to the target feature vector. A human-machine closed-loop benchmark path model is established using all historical input data within the historical benchmark time window and the corresponding actual target features. During parameter determination, the sum of the squares of the differences between the model output value and the actual target feature value is used as the first part of the computational cost, and the sum of the squares of all path parameter values ​​multiplied by 0.001 is used as the second part of the computational cost. The first part of the computational cost and the second part of the computational cost are added together to obtain the parameter solution objective. Each path parameter is adjusted until the parameter solution objective obtains the minimum value in this calculation process. For positions with a value of 0 in the causal reachability relation matrix, the corresponding path parameter remains 0 throughout the parameter solution process. Finally, all the obtained path parameters are arranged in a fixed matrix position to form the benchmark path parameter matrix. This baseline path parameter matrix is ​​used to describe the time-constrained path relationships between machine observations, human behavior, and machine control outputs under normal operating conditions.

[0018] In the above scheme, in S4, the current time window length is set to 2 seconds, and the window movement step is 100 milliseconds. Each time the time window moves, the prediction, residual calculation, and joint mutation judgment are re-performed.

[0019] For the first feature vector within the current time window, extract its historical features that satisfy causal reachability within the previous 500 milliseconds, arrange them in the same order as during the baseline path model training phase, and then calculate the predicted value of the first feature vector using the baseline path parameter matrix. Subtract the corresponding predicted value from the actual value of the first feature vector to obtain the first residual vector. Arrange all the first residual vectors in the current time window according to the sampling time to form a sequence of first residual vectors. For the second feature vector within the current time window, use the same process to calculate the predicted value of the second feature vector, subtract the corresponding predicted value from the actual value of the second feature vector to obtain the second residual vector, and form a sequence of second residual vectors according to the sampling time. Calculate the baseline mean and baseline covariance of the first residual using the first residual vectors in the historical baseline time window, and calculate the baseline mean and baseline covariance of the second residual using the second residual vectors in the historical baseline time window. For each set of first and second residuals in the current 2-second time window, calculate the standardized deviation of each residual relative to its own historical residual distribution. In the specific calculation, the current residual is first subtracted from the corresponding historical residual mean, and then scaled and normalized using the corresponding historical residual covariance to obtain the first and second residual deviation values. The first and second residual deviation values ​​at the same time point are added together, and then the arithmetic mean of all the sums within the current 2-second time window is calculated. This value is used as the joint mutation index for that time window. A historical joint mutation index sequence is generated using the same process for the historical baseline time window. The sum of all historical joint mutation indices is divided by the number of historical joint mutation indices to obtain the baseline mean of the joint mutation index; then the standard deviation of all historical joint mutation indices relative to this baseline mean is calculated. Adding three times the standard deviation to the baseline mean of the joint mutation index yields the event trigger threshold.

[0020] An adversarial event is identified when the joint mutation index of three consecutive current time windows reaches or exceeds the event trigger threshold, and the start time of the first time window of the three consecutive time windows is determined as the anchor point of the adversarial event. If the joint mutation index of any time window does not reach the event trigger threshold, no consecutive triggering result is generated.

[0021] In the above scheme, in S5, the adversarial event anchor point is used as the time boundary. The first residual vector sequence and the second residual vector sequence that are truncated for 1 second before the adversarial event anchor point are taken as the pre-event data, and the first residual vector sequence and the second residual vector sequence that are truncated for 1 second after the adversarial event anchor point are taken as the post-event data.

[0022] In the pre-event data, using the sampling time of each second residual vector as a matching benchmark, the first residual vector with the smallest absolute difference in sampling times is found from the first residual vectors. When the absolute difference in sampling times between the two is no more than 25 milliseconds, they are grouped into a cross-source residual sample; when it exceeds 25 milliseconds, the second residual vector is not used to establish a cross-source residual sample. The same time matching process is performed on the post-event data.

[0023] For each residual dimension in the first residual vector and each residual dimension in the second residual vector, calculate the linear correlation coefficient between them within 1 second before the event. The calculation process is as follows: calculate the arithmetic mean of the two residual dimensions within the window before the event; subtract the corresponding arithmetic mean from the residual value at each sampling point; multiply and sum the mean-free results of the two residual dimensions at the same sampling position; then divide the sum by the product of the square root of the sum of the squared mean-free results of the two residual dimensions. Arrange all correlation coefficients in fixed row and column positions to obtain the cross-source residual correlation matrix before the event.

[0024] The same calculation process is performed on the cross-source residual samples within 1 second after the event to obtain the post-event cross-source residual correlation matrix. Then, the value at each position in the post-event cross-source residual correlation matrix is ​​subtracted from the value at the same position in the pre-event cross-source residual correlation matrix to form the cross-source residual increment matrix. The absolute value of each value in the cross-source residual increment matrix is ​​taken. The confounding feature determination threshold is fixed at 0.30. When the absolute value of the cross-source residual increment corresponding to a certain position reaches or exceeds 0.30, the first and second residual dimensions corresponding to that position are determined as a post-event confounding feature pair; when the absolute value of the cross-source residual increment is less than 0.30, the corresponding residual dimension combination is not determined as a post-event confounding feature.

[0025] For each identified post-event confounding feature pair, its corresponding first and second residual components are read. Standardization is first performed using the arithmetic mean and standard deviation of the corresponding residual components within the historical reference time window; that is, the current residual component is subtracted from the historical arithmetic mean, and then divided by the historical standard deviation. Then, the standardized first residual component, the standardized second residual component, and the corresponding cross-source residual increment are multiplied to obtain the feature value of the post-event confounding feature at the current time. All post-event confounding features that meet the 0.30 criterion are arranged in order of first and second residual dimensions to form the confounding feature vector at the current time. All confounding feature vectors are arranged according to the sampling time to form a confounding feature vector sequence.

[0026] In the above scheme, in step S6, path parameters from human behavior to machine state and path parameters from machine control output to machine state are extracted from the baseline path parameter matrix. For the current moment, human behavior features that satisfy causal reachability within the previous 500 milliseconds are extracted. Each human behavior feature is multiplied by its corresponding human behavior path parameter, and the products acting on the same machine state dimension are accumulated to obtain the first state contribution component corresponding to the human behavior. Similarly, machine control output features that satisfy causal reachability within the previous 500 milliseconds are extracted. Each machine control output feature is multiplied by its corresponding machine control output path parameter, and the products acting on the same machine state dimension are accumulated to obtain the second state contribution component corresponding to the machine control output.

[0027] For each post-event confounding feature, first determine the human-side behavior residual dimension corresponding to the confounding feature, and then read the path parameters of the human-side behavior dimension acting on the machine state from the baseline path parameter matrix. Multiply the corresponding cross-source residual increment value by the above path parameters, and then divide by the sum of the square root of the sum of squares of the components of the above path parameters and one millionth to obtain the state mapping parameter corresponding to the confounding feature. Organize all state mapping parameters according to the order of the confounding features to form the mapping relationship from confounding features to machine states. Multiply each confounding feature value at the current time by the corresponding state mapping parameter, and accumulate the results acting on the same machine state dimension to obtain the third state contribution component corresponding to the post-event confounding factor. Then determine the attribution coefficients of the three types of paths: human-side behavior, machine control output, and post-event confounding factor. The values ​​of the three types of attribution coefficients are all limited to between 0 and 1. First, the baseline state component of the current machine state is calculated using the baseline path model. Then, the contribution components of the first state, the second state, and the third state are multiplied by their respective attribution coefficients. Finally, the baseline state component, the weighted contribution of human behavior, the weighted contribution of machine control output, and the weighted contribution of post-event confounding factors are added together to obtain the reconstructed machine state.

[0028] The actual machine state is subtracted from the reconstructed machine state dimension by dimension, and the difference in each dimension is squared. All squared results are then summed to obtain the sum of squares of the state reconstruction error at a single moment. The sum of squares of the state reconstruction errors at all moments within the current 2-second time window is accumulated, and then multiplied by 0.001 times the sum of squares of the three types of attribution coefficients. This result is used as the target for attribution coefficient calculation. The three types of attribution coefficients are adjusted while ensuring they are all within the range of 0 to 1, until the target for attribution coefficient calculation reaches its minimum value in this calculation process. This determines the attribution coefficients for human behavior, machine control output, and post-event confounding factors in the current time window. The state reconstruction error is calculated based on the Euclidean distance between the actual machine state and the reconstructed machine state, and then divided by the sum of the square root of the sum of squares of each dimension of the actual machine state and one millionth, thus obtaining the normalized state reconstruction error. This normalized state reconstruction error is used for subsequent path reliability calculations.

[0029] In the above scheme, in S7, the first state contribution component is multiplied by the human behavior attribution coefficient to obtain the corrected human behavior state contribution; the second state contribution component is multiplied by the machine control output attribution coefficient to obtain the corrected machine control output state contribution; and the third state contribution component is multiplied by the post-event confounding factor attribution coefficient to obtain the confounding state compensation amount. The generation process of the causal corrected state is as follows: first, the baseline state component output by the baseline path model is read; then, the corrected human behavior state contribution is added; then, the corrected machine control output state contribution is added; finally, the confounding state compensation amount is subtracted to obtain the causal corrected state at the current moment. All causal corrected states are arranged in chronological order to obtain the causal corrected state time series data. The causal corrected state time series data is input into the control optimization module. The control optimization module adopts model predictive control, predicting 20 machine control cycles each time, with each machine control cycle being 20 milliseconds, thus a single prediction covers 400 milliseconds. During the control optimization process, the sum of squares of the differences between the causal correction state and the target state in each dimension is calculated for each predictive control cycle. At the same time, 0.1 times the sum of squares of the values ​​of each dimension of the machine control output is calculated. The results of the above two parts are added together for 20 predictive control cycles to obtain the control cost.

[0030] The minimum and maximum permissible values ​​for each control output dimension are pre-stored at the machine control interface. The control optimization module adjusts the machine control outputs for 20 predictive control cycles, ensuring all control outputs fall between their corresponding minimum and maximum permissible values, to minimize the control cost during the current solution process. After each solution is completed, only the first control output of the 20 predictive control cycles is sent to the machine actuator. After 20 milliseconds, the machine observation data, human behavior data, and machine control output data are reacquired, and the state update and control optimization are performed again.

[0031] In the above scheme, in S8, once the adversarial event anchor point is determined in S4, the baseline path parameter matrix currently in use is immediately recorded, and the modification of the baseline path parameter matrix using data generated after the adversarial event anchor point is stopped, thus entering the parameter freeze state. In the parameter freeze state, data after the adversarial event anchor point continues to participate in the prediction of the first feature vector, the prediction of the second feature vector, residual calculation, cross-source residual correlation calculation, confounding feature extraction, state contribution calculation, path attribution, causal correction state generation, and machine control instruction calculation, but does not participate in the re-estimation and incremental update of the baseline path parameters.

[0032] During parameter freezing, the path reliability for each current time window is calculated based on the path attribution coefficient and state reconstruction error. First, the arithmetic mean of all normalized state reconstruction errors within the historical baseline time window is calculated as the historical reconstruction error baseline value. When calculating the current path reliability, the current normalized state reconstruction error is first divided by the sum of the historical reconstruction error baseline value and one part per million to obtain the first reliability decay. Then, the post-event confounding factor attribution coefficient is used as the second reliability decay. The first and second reliability decays are added together, and the negative of the sum is taken. The natural constant is then calculated using this negative number as the exponent, and the resulting value is used as the path reliability of the current data. The path reliability value is limited to between 0 and 1. Simultaneously, the recovery threshold is determined using the baseline mean and standard deviation of the historical joint mutation index. The recovery threshold is equal to the baseline mean of the historical joint mutation index plus 1.5 times the standard deviation of the historical joint mutation index.

[0033] During parameter freezing, a recovery status check is performed every 100 milliseconds. The parameter freeze is lifted only when the joint mutation index is less than the recovery threshold for five consecutive time windows, and the attribution coefficient of post-event confounding factors is less than 0.10 for five consecutive time windows. If either of the five consecutive windows fails to meet both conditions simultaneously, the consecutive count is reset to zero, and counting restarts from the next time window. After parameter freezing is lifted, the baseline path parameter matrix is ​​not directly updated using all new data. Path confidence is calculated for each new sample, with the path confidence update threshold fixed at 0.85. When a sample's path confidence reaches or exceeds 0.85, the sample is designated as a reliable update sample; when the path confidence is less than 0.85, the sample only participates in online control calculations and does not participate in baseline path parameter updates.

[0034] For reliable update samples, the baseline path parameter matrix is ​​incrementally updated using a recursive least squares approach with a forgetting factor of 0.995. Specifically, the parameter matrix and parameter covariance matrix saved after the last update are read; the historical input data of the current reliable update sample is multiplied by the parameter covariance matrix; then, the weighted energy corresponding to the historical input data under the current parameter covariance matrix and the forgetting factor of 0.995 are used to determine the parameter update gain for this update. Subsequently, the prediction result of the current reliable update sample is calculated using the baseline path parameter matrix before the update. The actual eigenvalues ​​are subtracted from the prediction result to obtain the current prediction error. This current prediction error is then multiplied by the parameter update gain, and the resulting parameter correction is added to the original baseline path parameter matrix. After the parameter matrix update is completed, the parameter covariance matrix is ​​corrected based on the same parameter update gain, and the covariance scaling is adjusted using a forgetting factor of 0.995. After each incremental update, all path parameters are checked again using the causal reachability matrix. Path parameters corresponding to positions with a value of 0 in the causal reachability matrix are directly set to 0, while positions with a value of 1 retain the result of the current incremental update. This ensures that the incremental learning process will not establish new paths that violate the constraints of sampling time, sampling period, and transmission delay.

[0035] Through the above processing, a human-machine closed-loop baseline path is established using 60 seconds of continuous normal operation data before the adversarial event occurs. When the adversarial event occurs, the joint mutation results of three consecutive current time windows are used to determine the adversarial event anchor point. After the anchor point is determined, post-event confounding features are extracted using the cross-source residual correlation changes one second before and after the anchor point. Three types of state contributions are calculated based on human behavior, machine control output, and post-event confounding factors, and the three types of path attribution coefficients are determined using the actual state reconstruction results. Subsequently, the confounding state compensation amount determined according to the attribution coefficients is deducted from the machine state to generate a causal correction state, and machine control commands are output accordingly. Simultaneously, the baseline path parameters are frozen from the adversarial event anchor point until the joint mutation index and the post-event confounding factor attribution coefficients simultaneously meet the recovery conditions for five consecutive time windows. After unfreezing, incremental updates are performed only using samples with a path confidence level of 0.85, thereby preventing post-event confounding data from being directly written into the human-machine closed-loop baseline path model.

[0036] The simulation platform of this invention: The experiment was completed in the MATLAB / Simulink (R2023b) desktop environment. The control algorithm ran with a fixed step size (20ms control cycle). The human side operation input was generated by the playback of the preset script plus the manually recorded joystick / pedal signal. The whole process does not require special hardware and can be reproduced on an ordinary workstation.

[0037] Controlled object and human-machine closed-loop structure: The controlled object is a two-dimensional planar trajectory tracking task (position x and position y) of a six-degree-of-freedom robotic arm end effector. The reference trajectory is a smooth sinusoidal curve with an amplitude of approximately 1m and a period of approximately 60s. The human side provides collaborative guidance through joystick displacement and direction inputs, and the machine control output is a driving force / speed command, forming a closed loop of human behavior—machine control output—machine observation state.

[0038] Data Acquisition and Parameter Setting: Three data sources: Machine observation data (position, velocity, acceleration, attitude, target distance, target relative velocity, sensor values, 12 channels in total) with a sampling period of 20ms; Human behavior data (joystick displacement, direction input, button, pedal, touch, 8 channels in total) with a sampling period of 50ms; Machine control output data with a sampling period of 20ms; Each data record includes the sampling time and reception time, with transmission delay randomly generated within 0–30ms, without interpolation or hard alignment, and a causal reachability matrix is ​​constructed according to S2 (time difference upper limit 500ms, lag order rounded up by the source sampling period); Feature extraction is performed by concatenating the current sampled value, adjacent rate of change, and the mean of 5 consecutive points, according to S1. The baseline path model was trained using normal data for 60 seconds after startup according to S3, with a regularization coefficient of 0.001. like Figure 6 As shown, the current time window is 2s, the step size is 100ms (S4); the event trigger threshold is the historical joint mutation mean + 3σ, the recovery threshold is the mean + 1.5σ, the anchoring is triggered for 3 consecutive windows, and the unfreezing condition is met for 5 consecutive windows; the confounding judgment threshold is 0.30 (S5); the attribution coefficient is limited to 0 and 1; the MPC prediction is 20 steps × 20ms (S7); the incremental update forgetting factor is 0.995, and the path confidence update threshold is 0.85 (S8).

[0039] Adversarial event injection method: Starting at t=60s and lasting for 2s, a band-limited random force perturbation (simulating sudden external adversarial / collision) is applied to the end effector of the robotic arm, causing an abnormal shift in the machine's observed state that is out of sync with the human input; the perturbation intensity is divided into three levels: low (amplitude 5N), medium (10N), and high (20N). During the perturbation, the human input remains unchanged in its normal collaborative guidance, thus forming a typical adversarial scenario where the normal correlation between human behavior and machine state is disrupted. Each experiment lasts for a total of 120s, and each working condition is independently repeated 5 times and the average is taken.

[0040] Comparison method settings: The method of this invention: complete execution of S1–S8 (causal correction state + freeze + trusted update); Freeze without correction: Perform event identification and model freezing, but control directly uses the uncorrected state (no S6–S7 attribution correction); Full data online updates: No freezing is performed; all data is continuously updated online after an event to update the baseline path model, and the causal correction state is controlled to not be used.

[0041] Experiment 1: Adversarial Event Recognition and Causal Correction Effects Experimental Procedure: After system startup, 60 seconds of normal data were collected to establish a baseline path model; a moderate-intensity disturbance pulse was injected at t=60 seconds; after t=60 seconds, the joint mutation index was calculated using the sliding window method S4, and anchoring was performed after three consecutive windows exceeded the threshold; the residuals before and after the anchor point were taken in step S5 to calculate the cross-source residual increment matrix and extract confounding features; the contributions and attribution coefficients of the three types of states were calculated in steps S6–S7, and the causal correction state was generated and sent to the MPC output control command; the model was frozen in step S8, and after the disturbance ended, it was unfrozen according to the recovery conditions and updated with a sample increment with a confidence level ≥0.85.

[0042] Table 1. Results of adversarial event identification and profanity feature extraction (medium intensity, 5 replicates)

[0043] The joint mutation index exceeded the threshold of 4.60 within one window after the start of the disturbance and was confirmed after three consecutive windows; the anchor point positioning deviation was less than 0.3s; the cross-source residual increment matrix (see...) Figure 6 Only the two injected and coupled cross-source channels exceeded the 0.30 threshold, while the increments of the remaining channels were within the normal fluctuation range of 0.05–0.14. This indicates that the correlation differences before and after the event can accurately and without redundancy locate the confounding features, verifying the effectiveness of S4–S5. The causal correction state reduced the peak deviation during the event period from 0.52m to 0.06m, verifying that path attribution correction can significantly reduce the probability of disturbances being misjudged as human intentions or normal machine responses.

[0044] Figure 3 This is the trajectory and joint mutation index curve of the present invention. Figure 3 The upper subplot shows the change of position dimension over time in a typical repeated experiment: the gray dashed line represents the reference trajectory, the red dotted line represents the uncorrected state, and the blue solid line represents the causal correction state of the method of this invention. Between t=60s (the adversarial event anchor point) and 62s (the end of the perturbation), the red curve deviates from the reference trajectory by a maximum of about 0.52m, while the blue curve deviates by only about 0.06m. Moreover, after the perturbation ends, the blue curve smoothly returns to the reference trajectory without significant overshoot. The lower subplot shows the joint mutation index: during the normal period, the index fluctuates steadily between 0 and 1. At the beginning of the perturbation, a peak of 8.6 appears, which is much higher than the event trigger threshold of 4.60. After the perturbation ends, the index quickly falls back and drops below the recovery threshold of 2.80, consistent with the unfreezing judgment process of S8.

[0045] Experiment 2: Comparison of tracking performance under different adversarial intensities Under three levels of perturbation intensity—low (5N), medium (10N), and high (20N)—the method of this invention, the freeze-only strategy, and the full data online update strategy were run, with each combination repeated 5 times. The RMSE of the end position tracking within 30 seconds (62–92 seconds) after the perturbation ended was statistically analyzed.

[0046] Table 2. Position tracking RMSE (cm, mean of 5 repetitions, values ​​in parentheses are minimum–maximum) under different adversarial intensities

[0047] At all three intensities, the RMSE of the method described in this invention was significantly lower than that of the two comparative methods. At medium intensity, it was 67.6% lower (7.1→2.3cm) than freezing without correction and 84.5% lower (14.8→2.3cm) than full data online update. As the intensity increased, the RMSE of the method described in this invention only increased slowly from 1.8cm to 3.0cm, while that of full data online update increased from 8.3cm to 24.1cm, with a degradation rate approximately 5 times that of the method described in this invention. This indicates that the causal correction state is the main source of the performance advantage, and the method has good robustness to disturbance intensity. Figure 4 In the dumbbell graph, the dots represent the mean and the horizontal bars represent the minimum-maximum range of 5 repetitions. It can be seen that the method of the present invention not only has a low mean, but also a significantly smaller spread of repeated experiments (error bar length), indicating that the control output is more stable and there are fewer human-machine conflicts.

[0048] Experiment 3: The impact of model freezing and reliable updates on the recovery process Experimental procedure: Under the same working conditions as Experiment 2, two indicators were recorded: (1) recovery time - the time required from the anchor point of the adversarial event until the position error is less than 5cm for 2 consecutive seconds and the joint mutation index is lower than the recovery threshold (2.80) for 5 consecutive windows; (2) steady-state RMSE 20s (92–112s) after recovery, which is used to check whether the data after the event has polluted the baseline path model.

[0049] Table 3 Recovery time after adversarial events under different update strategies (s, mean of 5 repetitions)

[0050] Table 4. Steady-state RMSE and parameter contamination test after recovery (medium intensity, mean of 5 replicates)

[0051] Under medium intensity, the recovery time of the method of this invention is 3.4s, which is 58.5% shorter than freezing without correction (8.2s) and 73.0% shorter than full data online update (12.6s). Under high intensity, the advantage is further expanded (5.6s vs. 21.3s). More importantly, the steady-state performance after recovery is significant: full data online update directly writes contaminated data after the event into the baseline path model, resulting in an RMSE of 4.6cm after recovery, which is nearly 3 times that of the normal baseline, indicating significant parameter contamination; the method of this invention, after freezing and updating with samples with a confidence level ≥0.85, achieves an RMSE of 1.7cm after recovery, which is basically the same as the normal baseline of 1.6cm, verifying the effect of S8 in suppressing parameter contamination. Figure 5 The bar chart visually shows the monotonic difference in recovery time among the three strategies at three intensity levels, with the method of this invention consistently showing the lowest difference.

[0052] Reducing hard alignment errors: Experiment 1 shows that under the conditions of no interpolation alignment and random transmission delay (0–30ms), the cross-source residual correlation before the event remains stable (increment 0.05–0.14), and the causal reachability matrix effectively isolates unreachable connections, verifying the role of asynchronous timing featureization and causal reachability modeling.

[0053] Adversarial event identification and confounding feature extraction: 100% identification success rate, anchor point deviation ≤0.3s, 2 pairs of confounding feature pairs were all hit with 0 false alarms, verifying the accuracy of joint mutation anchoring and extraction of pre- and post-event correlation differences. Improved control accuracy and stability: Peak deviation during the event period decreased from 0.52m to 0.06m; RMSE at three intensity levels was reduced by 67%–84% compared to the comparison method, and the spread in repeated experiments was smaller, indicating a significant reduction in human-machine control conflict and closed-loop oscillation. Suppressing parameter contamination and accelerating recovery: Recovery time was shortened by 58%–73%, and the steady-state RMSE after recovery remained at the normal baseline level (1.7cm vs. 4.6cm for full data online update), verifying the protective effect of the freeze-reliable update mechanism on the long-term health of the model.

[0054] Example 2: In this embodiment, S1 includes: Machine observation data, human behavior data, and machine control output data are collected separately, and associated with the corresponding original sampling time and data type identifier; Based on the original sampling time, the three types of data are sorted and their respective sampling periods are retained to construct the first, second and third time series sequences; The first, second, and third time series sequences are extracted according to a preset sliding window, and corresponding sample matrices are constructed respectively. Normalize and differ the sample matrices to extract state values, rates of change, and window statistical features, forming corresponding feature vectors; Arrange the first, second, and third feature vector sequences according to the sampling times corresponding to each feature vector.

[0055] In this embodiment, machine observation data, human behavior data, and machine control output data are obtained from their respective data acquisition interfaces. Machine observation data includes position, velocity, acceleration, attitude, target distance, and sensor detection values; human behavior data includes joystick displacement, pedal travel, direction input, button status, and touch input; machine control output data includes driving force commands, braking force commands, steering commands, torque commands, and actuator displacement commands. The sampling period for both machine observation data and machine control output data is set to 20 milliseconds, while the sampling period for human behavior data is set to 50 milliseconds. Each time the acquisition interface acquires a data record, it reads the timestamp generated by the acquisition device at the time of the sampling action and writes the timestamp into the original sampling time field of the corresponding sampling record. Simultaneously, it writes a data type identifier according to the data source, where the data type identifier for machine observation data is 01, for human behavior data it is 02, and for machine control output data it is 03. Therefore, each sampling record contains the original sampling time, data type identifier, and corresponding data value, and subsequent processing directly determines the time relationship based on the original sampling time.

[0056] The identified data is processed according to its data type. Sampling records with data type identifier 01 are sorted from earliest to latest based on their original sampling time, forming the first time series. Sampling records with data type identifier 02 are sorted according to the same rule, forming the second time series. Sampling records with data type identifier 03 are sorted according to the same rule, forming the third time series. Sorting only adjusts the storage order of data within the sequence; it does not change the original sampling time or sampling period in the sampling records. If duplicate records with identical original sampling times exist under the same data type, the record received first by the data acquisition interface is retained, and the remaining records with the same timestamp are deleted, ensuring that one sampling time in the same time series corresponds to only one set of sample values. After sorting, the time interval between adjacent sampling records in the first and third time series is 20 milliseconds, and the time interval between adjacent sampling records in the second time series is 50 milliseconds.

[0057] Three time series were each truncated using a sliding window containing five consecutive sampling points, with the window moving forward one sampling point at a time. For the first time series, one window contained five consecutive machine observation sampling points; for the second time series, one window contained five consecutive human behavior sampling points; and for the third time series, one window contained five consecutive machine control output sampling points. The five sampling points in the window were arranged in the order of their original sampling times, with each sampling point occupying one row in the sample matrix. The data channels were arranged in the columns of the sample matrix according to a pre-defined channel number order. Taking data containing seven machine observation channels as an example, the corresponding sample matrix consisted of 5 rows and 7 columns. The first row corresponded to the earliest sampling time within the window, the fifth row to the latest sampling time, and the first to seventh columns corresponded to the seven pre-determined machine observation channels. Human behavior data and machine control output data used the same matrix organization method. For each sample matrix, the original sampling time of the fifth row of sampling data was determined as the characteristic sampling time corresponding to that sample matrix.

[0058] After the sample matrix is ​​formed, it is normalized column by column according to the data channel. Each data channel has a pre-written lower and upper limit value consistent with the device interface range. For any data value in the sample matrix, the lower limit value of the corresponding channel is subtracted from the data value, and then the lower limit value is subtracted from the upper limit value of the corresponding channel to obtain the range difference. The difference between the data value and the lower limit value is then divided by the range difference, and the result is used as the normalized result for that data value. The same calculation process is used to process all elements in the sample matrix. Taking a velocity channel with a range of 0 to 40 meters per second and a current sample value of 10 meters per second as an example, the normalization process is as follows: first, 10 minus 0 equals 10; then, 40 minus 0 equals 40; finally, 10 divided by 40 equals 0.25. This converts data of different dimensions to a unified numerical scale while preserving the positional relationship of each sample value within its corresponding physical range.

[0059] Differential processing is performed based on the normalization result. For any data channel, the normalized value of the current sampling point is subtracted from the normalized value of the previous sampling point to obtain the change in adjacent sampling values. This change is then divided by the sampling time interval between the two sampling points to obtain the rate of change of the current sampling point. The divisor for machine observation data and machine control output data is 20 milliseconds, and the divisor for human behavior data is 50 milliseconds. For a window containing 5 sampling points, a total of 4 rates of change are calculated for the 2nd to 5th sampling points, and the rate of change of the 5th sampling point is used as the rate of change feature of the current window. The actual sampling time interval is used for division calculation, so that the rates of change formed by data from different sampling periods have corresponding time scales.

[0060] Window statistical features are calculated from five normalized sampled values ​​within the same sliding window. For each data channel, the five normalized values ​​of that channel within the current window are summed sequentially, and the sum is divided by 5 to obtain the window mean for that channel. Then, the differences between each of the five normalized values ​​and the window mean are calculated, the squares of each difference are summed, the sum is divided by 5, and the square root is calculated to obtain the window standard deviation for that channel. The window mean represents the state level of five consecutive sampled points, and the window standard deviation represents the dispersion of that data channel within the same window. The normalized value of the fifth sampled point in the current window is used as the state value, the rate of change corresponding to the fifth sampled point is used as the rate of change feature, and the window mean and window standard deviation are used as window statistical features. When forming the feature vector, the arrangement order is first determined according to the data channel number, and then the state value, rate of change, window mean, and window standard deviation are arranged sequentially within each data channel. The four types of features from all channels of the machine observation data are sequentially connected to form the first feature vector corresponding to the current window; the human behavior data forms the second feature vector in the same way; and the machine control output data forms the third feature vector in the same way. Assuming that the machine observation data contains 7 data channels, then one first feature vector contains 28 feature values; the feature vector dimensions of the human behavior data and the machine control output data are determined using the same calculation method according to the number of their respective data channels.

[0061] In this embodiment, for each feature vector generated, the original sampling time of the fifth sampling point in the sliding window that forms the feature vector is associated with that feature vector. Subsequently, the first, second, and third feature vectors are arranged from earliest to latest according to their original sampling times, forming a first feature vector sequence, a second feature vector sequence, and a third feature vector sequence. No time interpolation is performed on the three feature vector sequences, and different data sources are not resampled to the same time axis. Each feature vector retains its corresponding original sampling time, sampling period, and data type identifier. Through the above processing, the original sampled data is converted into state values, rates of change, and window statistical features with a definite time correspondence, providing direct input for subsequently determining the causal reachability relationship between feature vectors based on sampling time, sampling period, and transmission delay.

[0062] Example 3: Furthermore, S2 includes: Extract the sampling time, sampling period, and transmission delay corresponding to the first, second, and third time series to construct a time parameter vector; The causal reachability time interval is determined based on the sampling time difference, sampling period, and transmission delay between feature vectors from different data sources; The sampling time difference is matched with the causal reachability time interval. If the time reachability condition is not met, a zero-value reachability coefficient is generated. If the condition is met, a non-zero reachability weight is generated based on the deviation between the sampling time difference and the transmission delay. The reachability weights are arranged according to the indices of the source feature vector and the target feature vector to form a cross-source reachability submatrix; By combining the cross-source reachability sub-matrices according to the data source combination, the causal reachability matrix is ​​obtained.

[0063] In this embodiment, the first, second, and third feature vector sequences all retain their corresponding original sampling time, sampling period, and data reception time. For each feature vector generated, the data reception time corresponding to the final sampling point forming that feature vector is read, and the original sampling time is subtracted from the data reception time; the resulting time difference is determined as the transmission delay of that feature vector. The three time parameters are recorded in a fixed order of sampling time, sampling period, and transmission delay to form the time parameter vector corresponding to that feature vector. The first feature vector corresponds to machine observation data with a sampling period of 20 milliseconds; the second feature vector corresponds to human behavior data with a sampling period of 50 milliseconds; and the third feature vector corresponds to machine control output data with a sampling period of 20 milliseconds. Thus, each feature vector corresponds to a defined sampling time, sampling period, and transmission delay.

[0064] When determining the temporal reachability between different data sources, one feature vector is used as the source feature vector, and a feature vector from another data source is used as the target feature vector. First, the sampling time of the target feature vector is subtracted from the sampling time of the source feature vector to obtain the sampling time difference between the two. If the sampling time difference is less than 0, it indicates that the sampling of the target feature vector occurred before that of the source feature vector. This set of source and target feature vectors is not included in the subsequent reachability weight calculation, and its reachability coefficient is directly recorded as 0.

[0065] For a combination of feature vectors with a sampling time difference of not less than 0, the transmission delay of the source feature vector is used as the starting value of the causal reachability time interval. Then, the sampling period of the source feature vector is added to the sampling period of the target feature vector, and the sum is added to the aforementioned transmission delay to obtain the ending value of the causal reachability time interval. The starting value to the ending value constitutes the causal reachability time interval between this set of source and target feature vectors. This process uses the actual time required for data to be transmitted from the sampling end to the processing end as the earliest reachable time, while introducing the sampling periods of the two data sources into the termination boundary to absorb the time discrepancy caused by asynchronous sampling.

[0066] Taking the first feature vector pointing to the second feature vector as an example, the sampling period of the first feature vector is 20 milliseconds, and the sampling period of the second feature vector is 50 milliseconds. When the transmission delay of the corresponding source feature vector is 12 milliseconds, the starting value of the causal reachability time interval is 12 milliseconds. Adding 20 milliseconds and 50 milliseconds gives 70 milliseconds, and adding 70 milliseconds and 12 milliseconds gives 82 milliseconds. Therefore, the causal reachability time interval of this set of feature vectors is determined to be from 12 milliseconds to 82 milliseconds. After obtaining the causal reachability time interval, the sampling time difference is compared with the interval. When the sampling time difference is less than the starting value of the interval, it indicates that the source feature has not yet completed the prescribed time transmission process when the target feature is formed; when the sampling time difference is greater than the ending value of the interval, it indicates that the source feature has exceeded the effective time range corresponding to the target feature. In both of the above two cases, the corresponding reachability coefficient is recorded as 0. When the sampling time difference is between the starting value and the ending value of the interval and includes two boundary values, the set of feature vectors is considered to meet the time reachability condition, and a non-zero reachability weight is calculated.

[0067] The non-zero reachability weight is determined based on the deviation between the sampling time difference and the transmission delay of the source feature vector. The calculation involves first subtracting the transmission delay of the source feature vector from the sampling time difference and taking the absolute value of the difference to obtain the time deviation. Then, the sampling period of the source feature vector is added to the sampling period of the target feature vector to obtain a time scale reference value. Next, the time deviation is divided by the time scale reference value, and 1 is added to the result. Finally, 1 is divided by the result after adding 1, and the resulting value is the reachability weight of that set of feature vectors. According to this calculation method, when the sampling time difference equals the transmission delay, the time deviation is 0, and the reachability weight is 1. As the sampling time difference deviates from the transmission delay, the reachability weight decreases according to a defined calculation relationship and remains a non-zero value within the causal reachability time interval.

[0068] Specifically, based on the aforementioned 12-millisecond transmission delay, when the sampling time difference between the first and second feature vectors is 32 milliseconds, 32 milliseconds falls within the causal reachability time interval of 12 to 82 milliseconds. First, subtract 12 milliseconds from 32 milliseconds to obtain a 20-millisecond time deviation; then add the 20-millisecond sampling period of the first feature vector to the 50-millisecond sampling period of the second feature vector to obtain a 70-millisecond timescale baseline value; divide 20 by 70 to get 0.2857, add 1 to get 1.2857, and divide 1 by 1.2857 to obtain a reachability weight of 0.7778. If the sampling time difference is 90 milliseconds, this value exceeds the interval termination value of 82 milliseconds, and the corresponding reachability coefficient is directly recorded as 0. The above processing is performed on all feature vectors between the two data sources pairwise. When constructing the cross-source reachability submatrix, the index of the target feature vector is used as the matrix row index, and the index of the source feature vector is used as the matrix column index. The reachability weight calculated for each set of feature vectors is written into the corresponding row and column positions. Assuming the source data sequence contains 100 feature vectors and the target data sequence contains 60 feature vectors, a 60-row, 100-column cross-source reachability submatrix is ​​formed. A matrix element of 0 indicates that the corresponding source feature vector and target feature vector do not satisfy the time reachability condition; a matrix element with a non-zero value indicates that the time reachability condition is satisfied, and the non-zero value also represents the time reachability weight between the corresponding feature vectors.

[0069] For each directed data source combination, calculate six cross-source reachable sub-matrices corresponding to the following sequences: first eigenvector sequence to second eigenvector sequence, first eigenvector sequence to third eigenvector sequence, second eigenvector sequence to first eigenvector sequence, second eigenvector sequence to third eigenvector sequence, third eigenvector sequence to first eigenvector sequence, and third eigenvector sequence to second eigenvector sequence. Sub-matrices with different source data sources or different directions of action are calculated independently, without transposing or replacing sub-matrices with opposite directions. After completing the six cross-source reachable sub-matrices, determine the row block positions of the overall matrix based on the target data source and the column block positions based on the source data source. The first, second, and third data sources are respectively designated as row blocks 1, 2, and 3, and also as column blocks 1, 2, and 3, respectively. The six cross-source reachable sub-matrices are written into the off-diagonal block positions corresponding to their source and target data sources, with diagonal block positions corresponding to the same data source filled with 0. After concatenating all matrix blocks according to the above positions, a causal reachability matrix is ​​formed. When establishing the human-machine closed-loop baseline path model, only cross-source feature connections corresponding to non-zero elements in the causal reachability matrix are allowed to enter the path parameter estimation, and the value of the non-zero element is used as the time reachability weight of the corresponding connection.

[0070] Example 4: Furthermore, S3 includes: The first, second, and third feature vectors with reachability weights greater than zero in the causal reachability matrix are selected from the historical benchmark time window to form a benchmark feature sample set. The benchmark feature sample set is paired according to the reachability relationship between the source feature vector and the target feature vector to construct a cross-source feature sample matrix; Based on the cross-source feature sample matrix, a human-machine closed-loop baseline path model is established to represent the mapping relationship between human behavior, machine control output, and machine observation state. Using the reachable weights as sample weights, and aiming to minimize the weighted reconstruction error between the model-predicted features and the target features, the path parameters are solved. The path parameters are arranged according to the path direction between each data source to obtain the baseline path parameter matrix.

[0071] In this embodiment, the historical baseline time window is fixed at 60 seconds. When initially establishing the human-machine closed-loop baseline path model, the first, second, and third feature vectors are extracted from the first complete 60-second running interval after device initialization. After the model is put into online operation, the historical baseline time window is defined as the data from the 60 consecutive seconds prior to the current processing time that do not fall within the determined adversarial event interval. The obtained causal reachability matrix is ​​read, and the target feature vector is determined based on the row index, the source feature vector is determined based on the column index, and the corresponding reachability weights are read item by item. When a matrix element is equal to 0, the combination of the source and target feature vectors at that position is discarded; when a matrix element is greater than 0, the corresponding source and target feature vectors, their indices, and the value of the matrix element are saved. After completing the above filtering of all matrix elements within the historical baseline time window, a baseline feature sample set for baseline path modeling is obtained.

[0072] When performing cross-source pairing on the baseline feature sample set, the sample correspondence is directly established according to the source and target indices already determined in the causal reachability matrix, without re-performing temporal proximity matching. For each reachability weight greater than 0, the source feature vector corresponding to the column containing that weight is paired with the target feature vector corresponding to the row containing it to form a valid pairing sample, and the reachability weight is appended to this pairing sample. Subsequently, all valid pairing samples are arranged from earliest to latest according to the sampling time of the target feature vector; if the target sampling times are the same, they are arranged from smallest to largest according to the index of the source feature vector. One valid pairing sample occupies one row of the cross-source feature sample matrix. This row first writes all the eigenvalues ​​of the source feature vector, then writes all the eigenvalues ​​of the target feature vector, and records the reachability weight corresponding to this pairing sample at the end. The first, second, and third feature vectors form six cross-source feature sample matrices according to their direction of action, specifically including the first feature vector pointing to the second feature vector, the first feature vector pointing to the third feature vector, the second feature vector pointing to the first feature vector, the second feature vector pointing to the third feature vector, the third feature vector pointing to the first feature vector, and the third feature vector pointing to the second feature vector. Sample matrices are established for paths in opposite directions, without matrix transpose, so that each path retains its own corresponding sampling time difference, transmission delay, and reachability weight.

[0073] When establishing the human-machine closed-loop baseline path model, each path direction independently establishes a linear path mapping from source features to target features. For the first feature value in the target feature vector, each feature value of the source feature vector is multiplied by its corresponding path parameter, and then all products are summed. The resulting value is used as the predicted value of the target feature. The remaining features in the target feature vector are predicted item by item according to the same calculation process, ultimately forming a complete target predicted feature vector. The path from the second feature vector to the first feature vector represents the mapping of human behavior to the subsequent machine observation state, and the path from the third feature vector to the first feature vector represents the mapping of machine control output to the subsequent machine observation state. The remaining four path directions are used to describe the remaining cross-source mapping relationships between machine observation, human behavior, and machine control output. The six directed paths together constitute the human-machine closed-loop baseline path model.

[0074] When solving for path parameters, the determined reachability weights are directly used as the calculation weights for the corresponding paired samples. For any pair of paired samples, the target features of the sample are first predicted using the current path parameters. Then, the predicted feature value is subtracted from each feature value in the actual target features to obtain the dimension-wise reconstruction difference. The reconstruction differences are squared and summed to obtain the reconstruction error of the sample. The reconstruction error is then multiplied by the reachability weight of the sample to obtain the weighted reconstruction error of the sample. This calculation is performed on all valid paired samples in the same path direction, and all weighted reconstruction errors are summed. The result is used as the target parameter for solving the path direction.

[0075] In the specific solution process, the reachability weights of each sample group are first square-rooted. Then, the square root is multiplied by all source and target feature values ​​of that sample group. Following the sample arrangement, the source feature vectors processed by the square root of the weights are arranged row-by-row to form a weighted source feature matrix, and the target feature vectors processed by the same square root of the weights are arranged row-by-row to form a weighted target feature matrix. Subsequently, singular value decomposition is performed on the weighted source feature matrix. The reciprocal of each non-zero singular value is taken, and singular values ​​with a value of 0 are kept as 0. The matrix is ​​then recombinated according to the original positional relationships to obtain the generalized inverse matrix of the weighted source feature matrix. Multiplying this generalized inverse matrix by the weighted target feature matrix yields the path parameters for the current path direction. The path parameters obtained through this calculation process correspond to the minimum sum of weighted reconstruction errors for all valid samples.

[0076] Taking the path from the second feature vector to the first feature vector as an example, all human-side behavioral features with reachable weights greater than 0 within the historical reference time window are paired with machine observation features. Each group of human-side behavioral features is first multiplied by the square root of its corresponding reachable weight to form a weighted input sample; the paired machine observation features are multiplied by the square root of the same weight to form a weighted target sample. The weighted source feature matrix is ​​constructed from all weighted input samples, and the weighted target feature matrix is ​​constructed from all weighted target samples. Then, the path parameters from human-side behavior to machine observation state are obtained according to the aforementioned singular value decomposition and generalized inverse calculation process. The remaining 5 path directions are solved using the same parameter solution process, and the parameters are not shared between paths. After the 6 groups of path parameters are solved, the row blocks of the parameter matrix are determined according to the target data source, and the column blocks of the parameter matrix are determined according to the source data source. The first, second, and third data sources correspond to the 1st, 2nd, and 3rd row blocks, and also correspond to the 1st, 2nd, and 3rd column blocks, respectively. The 6 off-diagonal blocks between different data sources are written with the path parameters for the corresponding direction; the 3 diagonal blocks corresponding to the same data source are all filled with 0-value parameters. After concatenating the nine matrix blocks at fixed positions, the baseline path parameter matrix is ​​obtained. In subsequent predictions, the corresponding matrix blocks are directly read from the target and source data sources, and these blocks are used to perform mapping calculations on the current cross-source features, thereby obtaining the corresponding baseline path prediction results between the data sources.

[0077] Example 5: In this embodiment, S4 includes: Extract the second and third feature vectors that satisfy the causal reachability relationship within the current time window to form a first prediction input sequence, and extract the first and third feature vectors that satisfy the causal reachability relationship to form a second prediction input sequence; The first and second predicted input sequences are mapped using the reference path parameter matrix to obtain machine-side and human-side predicted feature vectors. The differences between the machine-predicted feature vector and the corresponding first feature vector, and between the human-predicted feature vector and the corresponding second feature vector are calculated respectively to obtain the first residual vector sequence and the second residual vector sequence; Based on the sliding window calculation of the residual amplitude changes and cross-source synchronization changes of the two residual vector sequences, a joint mutation index is generated. The sampling moment when the joint mutation index first reaches the preset mutation threshold is determined as the adversarial event anchor point.

[0078] In this embodiment, the current time window length is set to 2 seconds, and the window movement step size is set to 100 milliseconds. For each first feature vector within the current time window, the reachability weights of the second feature vector pointing to the first feature vector and the third feature vector pointing to the first feature vector are read according to the causal reachability matrix. Second and third feature vectors with reachability weights equal to 0 are not included in the current prediction process. Second and third feature vectors with reachability weights greater than 0 are arranged in order of sampling time, and their corresponding reachability weights are retained synchronously to form the prediction input corresponding to the first feature vector. The above processing is performed sequentially according to the sampling order of the first feature vectors within the current time window to form the first prediction input sequence.

[0079] For each second feature vector within the current time window, the reachable weights of the first feature vector pointing to that second feature vector and the third feature vector pointing to that second feature vector are read. Only the first and third feature vectors with reachable weights greater than 0 are retained and arranged according to the sampling time to form the prediction input corresponding to that second feature vector. After processing all second feature vectors within the current time window in sequence, the second prediction input sequence is obtained. Thus, the first prediction input sequence is used to predict the machine-side first feature vector, and the second prediction input sequence is used to predict the human-side second feature vector. When performing machine-side feature prediction, the path parameter blocks of the second feature vector pointing to the first feature vector and the path parameter blocks of the third feature vector pointing to the first feature vector are extracted from the baseline path parameter matrix, respectively. For a second feature vector in the first prediction input, each feature value in the feature vector is multiplied by the corresponding path parameter, and then the products are accumulated according to the target dimension of the first feature vector to obtain the machine-side path mapping result corresponding to the second feature vector; then the path mapping result is multiplied by the reachable weight corresponding to the second feature vector. The third feature vector is processed in the same way. The path parameter block pointing from the third feature vector to the first feature vector is used to obtain the machine-side path mapping result, and then multiplied by the corresponding reachability weight.

[0080] When a first feature vector to be predicted corresponds to multiple second and third feature vectors that satisfy causal reachability, all path mapping results corrected by reachability weights are accumulated according to the target feature dimension. Then, the accumulated result of each target dimension is divided by the sum of all reachability weights participating in this prediction to obtain the machine-side predicted feature vector at that moment. The calculation is performed sequentially according to the sampling time of the first feature vector to obtain the machine-side prediction result within the current time window.

[0081] The human-side prediction process maintains the same data processing relationship as the machine-side prediction process. Path parameter blocks pointing from the first feature vector to the second feature vector and from the third feature vector to the second feature vector are extracted from the baseline path parameter matrix. The first and third feature vectors in the second prediction input are multiplied and accumulated dimension-wise with their corresponding path parameters, and the resulting path mapping results are multiplied by their respective reachability weights. For the same second feature vector to be predicted, all weighted path mapping results are accumulated along the target dimension and then divided by the sum of the reachability weights participating in this prediction to obtain the human-side prediction feature vector at that sampling time. After obtaining the prediction result, the predicted value of the same dimension in the machine-side prediction feature vector is subtracted from each feature value in the actual first feature vector to obtain the first residual vector at that sampling time; all first residual vectors are arranged according to their sampling times to form a first residual vector sequence. The predicted value of the same dimension in the human-side prediction feature vector is subtracted from each feature value in the actual second feature vector to obtain the second residual vector, and a second residual vector sequence is formed according to the sampling times of the second feature vector.

[0082] In this embodiment, the residual amplitude is determined by the sum of the squares of all components in each residual vector. For a first residual vector, each residual component is squared, all squared results are summed, and the square root of the sum is taken as the machine-side residual amplitude at that sampling moment. The second residual vector is calculated using the same process to obtain the human-side residual amplitude. Taking the current sampling moment as the end, a 1-second window is extracted from the current window and the 1-second window immediately preceding the current window. The arithmetic mean of the machine-side residual amplitudes within the two windows is calculated. The arithmetic mean of the current 1-second window is subtracted from the arithmetic mean of the preceding 1-second window, and the absolute value of the difference is taken to obtain the change in the machine-side residual amplitude. The change in the human-side residual amplitude is calculated using the same process. Cross-source synchronization is determined based on the temporal correspondence between the first and second residual vector sequences. Using the sampling time of the second residual vector as a reference, the first residual vector with the smallest absolute difference in sampling times is selected from the first residual vector sequence. When the absolute difference in sampling times between the two residual vectors does not exceed 25 milliseconds, these two residual vectors are grouped into a cross-source synchronization sample; if the difference exceeds 25 milliseconds, the pairing is not established. After completing all residual pairings within the current 1-second window, the machine-side residual amplitude and human-side residual amplitude corresponding to each synchronization sample are read.

[0083] For the current 1-second window, first calculate the arithmetic mean of the machine-side residual amplitude and the human-side residual amplitude respectively. Then, subtract the machine-side arithmetic mean from each machine-side residual amplitude, and subtract the human-side arithmetic mean from the corresponding human-side residual amplitude. Multiply and sum the two mean-reduced values ​​for each group, and simultaneously calculate the sum of the squares of the mean-reduced machine-side and human-side values ​​respectively. Take the square root of each sum of squares and multiply them. Finally, divide the sum of the aforementioned multiplications by the product to obtain the cross-source synchronization coefficient for the current 1-second window. Perform the same calculation for the 1-second window immediately preceding the current window to obtain the cross-source synchronization coefficient for the previous 1-second window. Subtract the cross-source synchronization coefficient of the previous 1-second window from the current cross-source synchronization coefficient, and take the absolute value of the difference to obtain the cross-source synchronization change.

[0084] When generating the joint mutation index, the changes in machine-side residual amplitude, human-side residual amplitude, and cross-source synchronization amplitude are first obtained using the same calculation process within a historical benchmark time window. The arithmetic mean and standard deviation of these three types of historical data are then calculated for each. For the current time window, the standardized machine-side change value is obtained by subtracting its historical arithmetic mean from the current machine-side residual amplitude change value and then dividing by the corresponding historical standard deviation. The standardized human-side change value and the standardized synchronization change value are obtained using the same process for the human-side residual amplitude change and the cross-source synchronization change value, respectively. The absolute values ​​of the three standardized change values ​​are taken, summed, and divided by 3. The result is used as the joint mutation index for the current time window. The preset mutation threshold is determined based on the joint mutation index within the historical benchmark time window. The joint mutation index is calculated sequentially for each historical reference time window using a 100-millisecond moving step. The arithmetic mean of all historical joint mutation indices is obtained by summing them and dividing by the number of indices. The squares of the differences between each historical joint mutation index and this arithmetic mean are then calculated, summed, and divided by the number of historical joint mutation indices. Finally, the square root of the result is taken to obtain the standard deviation of the historical joint mutation index. The arithmetic mean of the historical joint mutation index is added to three times the standard deviation of the historical joint mutation index, and the resulting value is used as the preset mutation threshold. During online operation, the joint mutation index is calculated once every 100 milliseconds, and the sampling time corresponding to the end of the current 1-second sliding window is recorded as the joint mutation index. When the joint mutation index first reaches or exceeds the preset mutation threshold, the sampling time at the end of the window corresponding to that joint mutation index is determined as the adversarial event anchor point. Once the adversarial event anchor point is determined, it is not re-adjusted based on previous time positions; subsequent pre-event and post-event residual data use this sampling time as the time boundary.

[0085] Example 6: Specifically, S5 includes: Using the adversarial event anchor point as the boundary, the first residual vector sequence and the second residual vector sequence are respectively divided into the pre-event residual sequence and the post-event residual sequence; Based on the causal reachability matrix, the cross-source residual vectors are paired, and the outer product of the residual vectors is weighted and accumulated according to the corresponding reachability weights to construct the cross-source residual correlation matrix before the event and the cross-source residual correlation matrix after the event, respectively. The cross-source residual correlation matrix after the event is differiated from the cross-source residual correlation matrix before the event to obtain the cross-source residual increment matrix. Singular value decomposition is performed on the cross-source residual increment matrix, and the common residual direction whose corresponding singular values ​​reach the judgment threshold is determined as the post-event confounding feature; The post-event residual vector is projected onto the common residual direction, and the projection results are arranged according to the sampling time to form a mixed feature vector sequence.

[0086] The determination of the judgment threshold includes: Residual spectral features are extracted from multiple historical time windows corresponding to the cross-source residual correlation matrix before the event to form a historical singular value sequence. The pre-event residual strength benchmark is determined based on the preset quantile value of the historical singular value sequence; An adaptive judgment threshold is generated based on the pre-event residual strength benchmark and the preset safety margin; The singular values ​​corresponding to the current cross-source residual increment matrix are compared with the adaptive decision threshold, and the common residual directions that reach the adaptive decision threshold are determined as post-event confounding features.

[0087] In this embodiment, the aforementioned adversarial event anchor point is used as the time boundary for the residual data. The first and second residual vectors, occurring within one second before the adversarial event anchor point, are taken as the pre-event residual data. The first and second residual vectors, occurring within one second after the adversarial event anchor point, are taken as the post-event residual data. Both the first and second residual vectors retain the sampling time and index of the corresponding feature vector when the residual was generated, thus forming the pre-event first residual sequence, the pre-event second residual sequence, the post-event first residual sequence, and the post-event second residual sequence, respectively. When constructing the pre-event cross-source residual correlation matrix, the pre-event first and second residual vectors are paired according to the causal reachability matrix. For any one first residual vector and one second residual vector, the reachability weights in the causal reachability matrix are read according to the indices of their corresponding feature vectors. When the reachability weight is 0, no residual pairing is established; when the reachability weight is greater than 0, the first and second residual vectors are determined as a valid cross-source residual sample, and the corresponding reachability weights are recorded. The first residual vector and the second residual vector after the event are paired according to the same rules to obtain the effective cross-source residual samples after the event.

[0088] For each set of valid cross-source residual samples, calculate the outer product of the residual vectors. Specifically, multiply the first residual component of the first residual vector by all residual components of the second residual vector in sequence, and write the resulting products in the first row of the outer product matrix; then multiply the second residual component of the first residual vector by all residual components of the second residual vector in sequence, and write the products in the second row of the outer product matrix, until all components of the first residual vector have been calculated. The resulting outer product matrix has the same number of rows as the first residual vector and the same number of columns as the second residual vector.

[0089] Multiply all elements in each group of residual outer product matrices by the reachability weights corresponding to that group of residuals to obtain a weighted residual outer product matrix. For all valid cross-source residual samples before the event, sum the elements in the same row and column of each weighted residual outer product matrix; then calculate the sum of the reachability weights of all valid cross-source residual samples before the event, and divide the sum of the sums at each position by the sum of the reachability weights to obtain the cross-source residual correlation matrix before the event. The cross-source residual correlation matrix after the event is obtained using the same weighted summation and weight normalization process, so that the matrix before the event and the matrix after the event are on the same computational scale.

[0090] Subsequently, by subtracting the element at the same row and column position in the pre-event cross-source residual correlation matrix from each element in the post-event cross-source residual correlation matrix, and arranging the resulting differences in their original row and column order, the cross-source residual increment matrix is ​​obtained. Positive values ​​in the cross-source residual increment matrix correspond to an increase in the correlation of the cross-source residual combination after the event, while negative values ​​correspond to a reversal in the correlation direction. The absolute value of the value represents the change in the residual combination relative to the pre-event state.

[0091] Singular value decomposition (SVD) is performed on the cross-source residual increment matrix to extract common residual directions that appear frequently after the event. The calculation involves multiplying the cross-source residual increment matrix by its transpose, performing eigenvalue decomposition on the resulting symmetric matrix, and taking the square root of each non-negative eigenvalue to obtain the corresponding singular values, which are then arranged in descending order. The eigenvectors corresponding to each singular value serve as direction vectors in the first residual space. For each non-zero singular value, the transpose of the cross-source residual increment matrix is ​​multiplied by the corresponding first residual direction vector, and the result is divided by the singular value to obtain the corresponding direction vector in the second residual space. The corresponding first and second residual direction vectors together constitute a common residual direction. The threshold for determining whether the common residual direction belongs to post-event confounding features is determined based on historical residual data before the event. Twenty historical time windows, each 1 second long, are selected before the adversarial event anchor point, with a 100-millisecond interval between the start times of adjacent historical time windows. For each historical time window, a historical cross-source residual correlation matrix is ​​calculated according to the aforementioned process of causal reachability pairing, residual outer product, reachability weighting, and weight normalization. In chronological order, the cross-source residual correlation matrix corresponding to the previous historical time window is subtracted from the cross-source residual correlation matrix corresponding to the next historical time window to obtain the historical cross-source residual increment matrix. Then, the same singular value decomposition as the current data is performed on each historical cross-source residual increment matrix.

[0092] The singular values ​​obtained from each historical cross-source residual increment matrix are arranged in descending order of value and then grouped into historical singular value sequences according to the same order. The first-order historical singular value sequence consists of the first singular value in each historical time window, the second-order historical singular value sequence consists of the second singular value in each historical time window, and the remaining orders are formed in the same way. For any order of historical singular value sequence, the 20 singular values ​​are arranged in ascending order of value, and the 19th singular value after sorting is used as the pre-event residual strength benchmark for that order.

[0093] Specifically, a 20% safety margin is set based on the pre-event residual strength benchmark. This is achieved by multiplying the pre-event residual strength benchmark for each order by 0.20 to obtain the safety margin for that order. This safety margin is then added to the corresponding pre-event residual strength benchmark to obtain the adaptive decision threshold corresponding to the common residual direction for that order. For example, if the pre-event residual strength benchmark for a certain order is 0.50, the safety margin is 0.10, and the adaptive decision threshold for that order is 0.60. The decision threshold for each order of singular values ​​is calculated using a separate historical singular value sequence for that order; no single threshold is shared.

[0094] After the current cross-source residual increment matrix completes singular value decomposition, the current first-order singular value is compared with the first-order adaptive decision threshold, the current second-order singular value is compared with the second-order adaptive decision threshold, and so on for other orders. When the current singular value reaches or exceeds the adaptive decision threshold of the corresponding order, the first and second residual direction vectors corresponding to the singular value are retained, and the common residual direction formed by the two is determined as the post-event confounding feature direction; when the current singular value is less than the corresponding decision threshold, the common residual direction is discarded. When there are two or more common residual directions that meet the decision condition, they are arranged in descending order according to the corresponding current singular values. After determining the post-event confounding feature direction, the valid post-event cross-source residual samples are projected onto each common residual direction. For a common residual direction, each residual component of the first residual vector after the event is multiplied by the component of the corresponding first residual direction vector, and all products are summed to obtain the machine-side projection value. The paired second residual vector after the event is then multiplied and summed with the second residual direction vector in the same way to obtain the human-side projection value. Each common residual direction that meets the judgment criteria yields one machine-side projection value and one human-side projection value.

[0095] In this embodiment, the machine-side and human-side projection values ​​are arranged in descending order of the singular values ​​corresponding to the common residual direction to form the hybrid feature vector corresponding to the current cross-source residual sample. When there are two or more valid cross-source residual pairs at the same sampling time, the projection values ​​at the same feature position are multiplied by their respective reachable weights and summed. The summation result is then divided by the sum of the reachable weights of all valid pairs at that sampling time, thus obtaining a unique hybrid feature vector for that sampling time. Each hybrid feature vector uses the later sampling time in the corresponding cross-source residual sample as a time marker, and they are arranged in order of time markers to form a hybrid feature vector sequence. Therefore, the hybrid feature vector sequence only retains the common residual components that meet the adaptive judgment conditions relative to the historical residual spectrum before the event, and uses these as input for calculating the state contribution of subsequent hybrid factors after the event.

[0096] Example 7: Furthermore, S6 includes: The human-side induced control component that can be explained by the second eigenvector in the third eigenvector is calculated using the reference path parameter matrix, and the human-side induced control component is subtracted from the third eigenvector to obtain the machine autonomous control characteristics; Based on the baseline path parameter matrix, path mapping is performed on the second feature vector and the machine autonomous control feature to obtain the first state contribution component and the second state contribution component. Based on the mapping relationship between the hybrid feature vector sequence and the machine observation state, the third state contribution component is obtained, and the state contribution feature matrix is ​​constructed from the first, second and third state contribution components. The state reconstruction error is calculated based on the state contribution feature matrix and the actual machine state vector, and a constrained optimization model is established with the goal of minimizing the state reconstruction error. Solve the constrained optimization model under the condition that the path attribution coefficients are non-negative and satisfy normalization to obtain the path attribution coefficient vector.

[0097] In this embodiment, the machine control output in the third feature vector is first distinguished into control components caused by human behavior and control components generated autonomously by the machine. For the current third feature vector, the second feature vector with a reachability weight greater than 0 before its sampling time is found according to the causal reachability matrix, and the path parameters pointing from the second feature vector to the third feature vector in the reference path parameter matrix are read. For each second feature vector that meets the conditions, its feature values ​​are multiplied by the corresponding path parameters, and the results are accumulated according to the target dimensions of the third feature vector to obtain the control mapping result corresponding to the second feature vector; then the control mapping result is multiplied by the reachability weight corresponding to the second feature vector. When there are more than two second feature vectors, the weighted results of the same target dimension are added together and divided by the sum of all reachability weights involved in the calculation. The result is used as the human-induced control component in the current third feature vector. After obtaining the human-induced control component, the value at the same position in the human-induced control component is subtracted from each feature value in the third feature vector, and the difference is arranged according to the original feature order of the third feature vector to form the machine autonomous control feature. This process only deducts the control portion that can be explained by the baseline path from the second feature vector to the third feature vector, so that the remaining control features can be used to characterize the state action of the machine control side independent of the current human-side induced path.

[0098] The first state contribution component is obtained by mapping the second feature vector to the machine observation state. For the current first feature vector, the correspondence between the second feature vector and the first feature vector in the causal reachability matrix, where the reachability weight is greater than 0, is read. The path parameters between the second feature vector and the first feature vector in the baseline path parameter matrix are extracted. Each eigenvalue of each second feature vector is multiplied by its corresponding path parameter, summed according to the machine observation target dimension, and then multiplied by the corresponding reachability weight. When multiple valid second feature vectors exist, the weighted results for each target dimension are summed separately and then divided by the sum of the reachability weights used in this calculation. Subsequently, only the target dimension corresponding to the machine observation state value is extracted and arranged according to the machine state channel number to obtain the first state contribution component.

[0099] The second state contribution component is calculated using the aforementioned machine autonomous control features. The path parameters pointing from the third feature vector to the first feature vector are read from the baseline path parameter matrix. Each feature value in the machine autonomous control features is multiplied by the corresponding path parameter, and then accumulated according to the target dimension of the machine's observed state. The resulting state dimension values ​​are arranged in the same machine state channel order as the first state contribution component, forming the second state contribution component. Since the third feature entering this path mapping has already deducted the human-induced control component, the indirect effect formed by human behavior through machine control output will not be repeatedly counted as machine autonomous control action.

[0100] The third state contribution component establishes a mapping relationship between the corresponding samples of the promiscuous feature vector sequence and the machine observation state. Using the time stamp of each promiscuous feature vector as a reference, it searches for the first residual vector in the first residual vector sequence whose sampling time matches the first residual vector. If no identical sampling time exists, it selects the first residual vector with the smallest absolute time difference, not exceeding 20 milliseconds. The promiscuous feature vector is used as the input sample, and the residual component in the first residual vector corresponding to the machine state value is used as the target sample. The promiscuous mapping sample set is composed of all valid samples within one second after the event.

[0101] When obtaining the mapping parameters from profanity features to machine observation states, the profanity feature vectors in the profanity mapping sample set are arranged row by row according to the sampling time to form the input sample matrix; the corresponding machine state residuals are arranged row by row according to the same sample order to form the target sample matrix. Singular value decomposition is performed on the input sample matrix, taking the reciprocal of each non-zero singular value, keeping the singular values ​​of 0 as 0, and then combining them according to the original decomposition positions to obtain the generalized inverse of the input sample matrix. This generalized inverse is multiplied by the target sample matrix to obtain the mapping parameters from profanity features to machine observation states. These parameters minimize the sum of squares of the reconstructed machine state residuals for all samples within one second after the event. When mapping the current profanity feature vector to a state, each profanity feature value is multiplied by its corresponding mapping parameter, and the results are accumulated according to the machine state dimension. The accumulated results for each state dimension constitute the third state contribution component. The first, second, and third state contribution components are all arranged according to the same machine state channel, thus enabling direct subsequent state reconstruction.

[0102] Specifically, when constructing the state contribution feature matrix, the machine state channels are used as matrix rows, and the three types of state contributions are used as matrix columns. The first column contains the first state contribution component, the second column contains the second state contribution component, and the third column contains the third state contribution component. The actual machine state vector is extracted from the current first feature vector, selecting only the features corresponding to the current state value, excluding the rate of change and window statistics features, and arranged according to the row order of the state contribution feature matrix.

[0103] During state reconstruction, one path attribution coefficient is assigned to each of the three types of state contributions. The first state contribution component is multiplied by the first attribution coefficient, the second state contribution component by the second attribution coefficient, and the third state contribution component by the third attribution coefficient. Then, for each machine state dimension, the three weighted contribution values ​​are summed to obtain the corresponding reconstructed state value. The reconstructed state value is subtracted from the actual machine state value, and the difference for each state dimension is squared. The sum of all squared values ​​is divided by the machine state dimension, and the result is taken as the state reconstruction error. The constrained optimization model aims to minimize the state reconstruction error, and stipulates that all three path attribution coefficients are not less than 0, and the sum of the three path attribution coefficients is equal to 1. Since only three attribution coefficients are involved, this embodiment uses an effective combination enumeration method to solve the problem. Seven combinations are calculated sequentially: all three attribution coefficients participating simultaneously, any two attribution coefficients participating, and a single attribution coefficient participating. For any combination, the attribution coefficients of the unparticipated factors are fixed at 0. The inner product operation is performed on the state contribution column and the actual machine state vector, respectively. The linear stationary point equation of the state reconstruction error is established from the inner product result, and the constraint condition that the sum of the participating attribution coefficients is equal to 1 is added. The attribution coefficients corresponding to the combination are obtained by Gaussian elimination.

[0104] In this process, for each set of candidate attribution coefficients, it is first checked whether each coefficient is not less than 0; if there is a negative value, the result of that set is discarded. When the non-negativity condition is met, the contributions of the three types of states are re-weighted using the attribution coefficients of that set, and the corresponding state reconstruction error is calculated and recorded. After all seven combinations have been calculated, the set with the smallest state reconstruction error is selected from the candidate results that meet the non-negativity and normalization constraints, and arranged in the order of human behavior, machine autonomous control, and post-event confounding factors to obtain the path attribution coefficient vector. Through calculation, the human-induced control component in the third feature vector is first extracted separately and deducted from the machine control feature. Then, the three types of state contributions corresponding to human behavior, machine autonomous control, and post-event confounding factors are obtained respectively. The reconstruction error of the actual machine state is used as a unified constraint to determine the attribution ratio of the three types of paths, providing input for subsequent correction of the machine state according to the path attribution results.

[0105] Example 8: S7 includes: A path correction weight matrix is ​​constructed based on the path attribution coefficient vector, and the contribution components of the first, second and third states are corrected respectively. The remaining state components are calculated based on the corrected three state contribution components and the actual machine state vector, and the residual confidence is determined based on the amplitude stability of the remaining state components and their correlation with the hybrid feature vector. The remaining state components whose residual confidence reaches the confidence threshold are fused with the corrected first, second and third state contribution components to obtain the causal corrected state vector, and arranged according to the sampling time to form the causal corrected state time series data. The causal correction state timing data is input into the control optimization module, and the control quantity is calculated and machine control instructions are generated in combination with the current control constraints. Machine control commands are sent to the machine actuators to perform human-machine collaboration and adversarial control.

[0106] In this embodiment, path correction is based on the path attribution coefficient vector obtained above. The path attribution coefficient vector includes, in sequence, human behavior attribution coefficient, machine autonomous control attribution coefficient, and post-event confounding factor attribution coefficient. All three coefficients are not less than 0, and the sum of the three coefficients is 1. A 3x3 path correction weight matrix is ​​constructed, with the human behavior attribution coefficient written in the first row and first column, the machine autonomous control attribution coefficient in the second row and second column, and the post-event confounding factor attribution coefficient in the third row and third column. All other positions are assigned a value of 0. Each state value in the first state contribution component is multiplied by the human behavior attribution coefficient to obtain the corrected first state contribution component; each state value in the second state contribution component is multiplied by the machine autonomous control attribution coefficient to obtain the corrected second state contribution component; and each state value in the third state contribution component is multiplied by the post-event confounding factor attribution coefficient to obtain the corrected third state contribution component. The three correction results are arranged in the same machine state channel order, so that the same row position corresponds to the same machine state quantity.

[0107] The residual state component is determined based on the difference between the actual machine state and the contributions of the three types of attributed states. For a single machine state channel, the corrected contributions of the first, second, and third states are first added together. Then, the actual machine state value of that channel is subtracted from the summed value, and the difference is taken as the residual state value of that channel. All residual state values ​​are calculated and arranged sequentially according to the machine state channel number to form the residual state component at the current sampling time. This residual state component represents the state component that remains in the actual machine state after interpretation by the three types of identified paths. The residual confidence is determined based on the amplitude stability of the residual state component and its correlation with the profane feature vector. A continuous 1-second data segment is taken with the current sampling time as the end, and a window movement step of 100 milliseconds is used. For each residual state component within the window, each state value is squared, and the square root of the sum of all squared values ​​is taken to obtain the residual state amplitude at the corresponding sampling time. Add up all remaining state amplitudes within the current 1-second window and divide by the number of samples to obtain the average value of the remaining state amplitudes; then subtract this average value from each remaining state amplitude, square the difference, and sum them up. Divide the sum by the number of samples and take the square root to obtain the standard deviation of the remaining state amplitudes.

[0108] Furthermore, amplitude stability is calculated as the ratio of the standard deviation of the residual state amplitude to the average of the residual state amplitudes. First, the standard deviation of the residual state amplitude is divided by the sum of the average of the residual state amplitudes and one part per million. Then, 1 is added to the result, and finally, 1 is divided by the result to obtain the amplitude stability. When the residual state amplitudes are completely consistent within the current window, the standard deviation is 0, corresponding to an amplitude stability of 1. As the amplitude dispersion within the window increases, the amplitude stability decreases according to the above calculation relationship. The correlation between the residual state components and the profane feature vectors is calculated using the amplitude sequence within the same time window. First, the eigenvalues ​​of each profane feature vector are squared and summed. Then, the square root of the sum is taken to obtain the profane feature amplitude at that sampling time. Using the time stamp of the profane feature vector as a reference, the sample with the smallest absolute time difference (not exceeding 20 milliseconds) in the residual state amplitude sequence is selected for pairing. For all valid pairs within the current 1-second window, calculate the average value of the remaining state amplitude and the average value of the hybrid feature amplitude respectively. Then, subtract the average value from each remaining state amplitude and its corresponding hybrid feature amplitude.

[0109] The process involves multiplying the two mean-reduction results for each pair and summing all products. Simultaneously, the sum of the squares of the mean-reduction results for all remaining states and the sum of the squares of the mean-reduction results for all confounding features are calculated separately. The square roots of these two sums are multiplied, and the sum of the aforementioned products is divided by this product. The absolute value of the result is then taken to obtain the correlation degree, which is limited to 0 to 1. Subtracting the correlation degree from 1 yields the independence degree between the remaining state and the confounding factors. The amplitude stability and the independence degree are multiplied by 0.5, and the two results are summed to obtain the residual confidence level of the current remaining state component. Thus, the residual confidence level is constrained by both the amplitude fluctuation of the remaining state itself and its synchronization with the post-event confounding features. The confidence threshold is determined using 20 historical windows prior to the adversarial event anchor point, with each window being 1 second long and the start time interval between adjacent windows being 100 milliseconds. For each of the 20 historical windows, the residual confidence level is calculated using the process described above. The 20 residual confidence levels are summed and divided by 20 to obtain the average historical residual confidence level. The difference between each historical residual confidence level and this average is calculated, and the squared differences are summed. The sum is then divided by 20 and the square root is taken to obtain the standard deviation of the historical residual confidence level. The confidence threshold is obtained by subtracting twice the standard deviation from the average historical residual confidence level. If the result is less than 0, the confidence threshold is set to 0.

[0110] When the current residual confidence level reaches or exceeds the confidence level threshold, the current remaining state component is retained; when the current residual confidence level is less than the confidence level threshold, all state values ​​in the current remaining state component are set to 0. Thus, only remaining states that simultaneously meet the requirements of amplitude stability and confounding independence are included in the causal correction state calculation. When generating the causal correction state vector, the corrected first state contribution component, the corrected second state contribution component, and the remaining state component retained after confidence level determination are added item by item according to the machine state channel. The corrected third state contribution component is not added to this addition result as a valid state contribution; it has already been deducted from the actual machine state as a confounding path-corresponding state component during the remaining state calculation stage. The resulting state vector retains the human-side behavioral effects confirmed by path attribution, the machine's autonomous control effects, and the remaining states determined by confidence level determination, while excluding state contributions corresponding to confounding factors after the event.

[0111] For each causal correction state vector obtained, it is assigned the same sampling time as the current actual machine state and arranged in chronological order according to the sampling time to form causal correction state time-series data. The control optimization module performs a control calculation every 20 milliseconds and reads the causal correction states corresponding to the 20 consecutive control cycles before the current time. Therefore, each control calculation uses 400 milliseconds of continuous state data. The control optimization module adopts a prediction length of 20 control cycles. First, based on the current causal correction state, the target machine state, and the pre-calibrated machine state transition parameters, the predicted state corresponding to the next 20 control cycles is recursively calculated. For each prediction cycle, the difference between the predicted state and the target state in each state channel is calculated, and the squared differences are summed to obtain the state tracking error for that cycle. The control quantities of each control channel in that cycle are squared and summed, then multiplied by 0.1 to obtain the control quantity cost. The control quantity of the current prediction cycle is subtracted from the control quantity of the previous prediction cycle, the difference in each control channel is squared and summed, then multiplied by 0.05 to obtain the control change cost. The state tracking error, control cost, and control change cost of the 20 prediction cycles are accumulated and used as the target value for solving the control variable in this test.

[0112] The current control constraints are read from the machine control interface parameter table. For each control channel, the control quantity must not be less than the minimum allowable value recorded in the interface parameter table, nor exceed the maximum allowable value; the absolute value of the change in control quantity between two adjacent control cycles must not exceed the single-cycle allowable change of the corresponding channel; the predicted machine state values ​​must not exceed the upper and lower limits of state safety stored in the equipment parameter table. Under the above constraints, quadratic programming is used to solve for the control quantities corresponding to 20 prediction cycles, so that the aforementioned control target value is minimized. After this optimization is completed, the control quantity of the first control cycle is extracted from the control results of the 20 prediction cycles. The instruction fields corresponding to the driving force, braking force, steering amount, torque, or actuator displacement are written according to the channel order specified by the machine execution interface, and the current control cycle number and sending time are written at the same time to generate machine control instructions. The machine control instructions are sent to the corresponding actuators via the control bus, and the actuators complete the actions within the current 20-millisecond control cycle according to the received control values. After entering the next control cycle, the machine observation data, human behavior data, and machine control output data are reacquired, and the control instructions are recalculated based on the updated causal correction state, thereby continuously executing human-machine collaborative and adversarial control.

[0113] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A human-computer collaboration and adversarial control optimization method for complex dynamic environments, characterized in that, include: S1. Acquire machine observation data, human behavior data, and machine control output data; construct the first, second, and third time series according to the sampling time; and extract the first, second, and third feature vector sequences. S2. Determine the causal reachability relationship between feature vectors based on the sampling time, sampling period and transmission delay of each time series, and construct the causal reachability relationship matrix; S3. Based on the first, second and third eigenvectors that satisfy the causal reachability relationship within the historical benchmark time window, establish a human-machine closed-loop benchmark path model and obtain the benchmark path parameter matrix. S4. Predict the first and second eigenvectors within the current time window using the baseline path parameter matrix to obtain the first and second residual vector sequences, and determine the adversarial event anchor point based on their joint mutation. S5. Based on the adversarial event anchor points, construct the pre-event cross-source residual correlation matrix and the post-event cross-source residual correlation matrix respectively. Perform difference operation on the two to obtain the cross-source residual increment matrix, and extract the post-event hybrid features that reach the judgment threshold from it to form a hybrid feature vector sequence. S6. Based on the baseline path parameter matrix and the sequence of confounding feature vectors, determine the first, second and third state contribution components corresponding to human behavior, machine control output and post-event confounding factors, and determine the path attribution coefficient vector with state reconstruction error as a constraint. S7. Correct the contribution components of the first, second and third states according to the path attribution coefficient vector, generate causal correction state time series data, and input it into the control optimization module to generate machine control instructions to execute human-machine cooperation and adversarial control.

2. The human-machine collaboration and adversarial control optimization method for complex dynamic environments according to claim 1, characterized in that, S1 includes: Machine observation data, human behavior data, and machine control output data are collected separately, and associated with the corresponding original sampling time and data type identifier; Based on the original sampling time, the three types of data are sorted and their respective sampling periods are retained to construct the first, second and third time series sequences; The first, second, and third time series sequences are extracted according to a preset sliding window, and corresponding sample matrices are constructed respectively. Normalize and differ the sample matrices to extract state values, rates of change, and window statistical features, forming corresponding feature vectors; Arrange the first, second, and third feature vector sequences according to the sampling times corresponding to each feature vector.

3. The human-machine collaboration and adversarial control optimization method for complex dynamic environments according to claim 1, characterized in that, S2 includes: Extract the sampling time, sampling period, and transmission delay corresponding to the first, second, and third time series to construct a time parameter vector; The causal reachability time interval is determined based on the sampling time difference, sampling period, and transmission delay between feature vectors from different data sources; The sampling time difference is matched with the causal reachability time interval. If the time reachability condition is not met, a zero-value reachability coefficient is generated. If the condition is met, a non-zero reachability weight is generated based on the deviation between the sampling time difference and the transmission delay. The reachability weights are arranged according to the indices of the source feature vector and the target feature vector to form a cross-source reachability submatrix; By combining the cross-source reachability sub-matrices according to the data source combination, the causal reachability matrix is ​​obtained.

4. The human-machine collaboration and adversarial control optimization method for complex dynamic environments according to claim 1, characterized in that, S3 includes: The first, second, and third feature vectors with reachability weights greater than zero in the causal reachability matrix are selected from the historical benchmark time window to form a benchmark feature sample set. The benchmark feature sample set is paired according to the reachability relationship between the source feature vector and the target feature vector to construct a cross-source feature sample matrix; Based on the cross-source feature sample matrix, a human-machine closed-loop baseline path model is established to represent the mapping relationship between human behavior, machine control output, and machine observation state. Using the reachable weights as sample weights, and aiming to minimize the weighted reconstruction error between the model-predicted features and the target features, the path parameters are solved. The path parameters are arranged according to the path direction between each data source to obtain the baseline path parameter matrix.

5. The human-machine collaboration and adversarial control optimization method for complex dynamic environments according to claim 1, characterized in that, S4 includes: Extract the second and third feature vectors that satisfy the causal reachability relationship within the current time window to form a first prediction input sequence, and extract the first and third feature vectors that satisfy the causal reachability relationship to form a second prediction input sequence; The first and second predicted input sequences are mapped using the reference path parameter matrix to obtain machine-side and human-side predicted feature vectors. The differences between the machine-predicted feature vector and the corresponding first feature vector, and between the human-predicted feature vector and the corresponding second feature vector are calculated respectively to obtain the first residual vector sequence and the second residual vector sequence; Based on the sliding window calculation of the residual amplitude changes and cross-source synchronization changes of the two residual vector sequences, a joint mutation index is generated. The sampling moment when the joint mutation index first reaches the preset mutation threshold is determined as the adversarial event anchor point.

6. The human-machine collaboration and adversarial control optimization method for complex dynamic environments according to claim 1, characterized in that, S5 includes: Using the adversarial event anchor point as the boundary, the first residual vector sequence and the second residual vector sequence are respectively divided into the pre-event residual sequence and the post-event residual sequence; Based on the causal reachability matrix, the cross-source residual vectors are paired, and the outer product of the residual vectors is weighted and accumulated according to the corresponding reachability weights to construct the cross-source residual correlation matrix before the event and the cross-source residual correlation matrix after the event, respectively. The cross-source residual correlation matrix after the event is differiated from the cross-source residual correlation matrix before the event to obtain the cross-source residual increment matrix. Singular value decomposition is performed on the cross-source residual increment matrix, and the common residual direction whose corresponding singular values ​​reach the judgment threshold is determined as the post-event confounding feature; The post-event residual vector is projected onto the common residual direction, and the projection results are arranged according to the sampling time to form a mixed feature vector sequence.

7. The human-machine collaboration and adversarial control optimization method for complex dynamic environments according to claim 1, characterized in that, S6 includes: The human-side induced control component that can be explained by the second eigenvector in the third eigenvector is calculated using the reference path parameter matrix, and the human-side induced control component is subtracted from the third eigenvector to obtain the machine autonomous control characteristics; Based on the baseline path parameter matrix, path mapping is performed on the second feature vector and the machine autonomous control feature to obtain the first state contribution component and the second state contribution component. Based on the mapping relationship between the hybrid feature vector sequence and the machine observation state, the third state contribution component is obtained, and the state contribution feature matrix is ​​constructed from the first, second and third state contribution components. The state reconstruction error is calculated based on the state contribution feature matrix and the actual machine state vector, and a constrained optimization model is established with the goal of minimizing the state reconstruction error. Solve the constrained optimization model under the condition that the path attribution coefficients are non-negative and satisfy normalization to obtain the path attribution coefficient vector.

8. The human-machine collaboration and adversarial control optimization method for complex dynamic environments according to claim 1, characterized in that, S7 includes: A path correction weight matrix is ​​constructed based on the path attribution coefficient vector, and the contribution components of the first, second and third states are corrected respectively. The remaining state components are calculated based on the corrected three state contribution components and the actual machine state vector, and the residual confidence is determined based on the amplitude stability of the remaining state components and their correlation with the hybrid feature vector. The remaining state components whose residual confidence reaches the confidence threshold are fused with the corrected first, second and third state contribution components to obtain the causal corrected state vector, and arranged according to the sampling time to form the causal corrected state time series data. The causal correction state timing data is input into the control optimization module, and the control quantity is calculated and machine control instructions are generated in combination with the current control constraints. Machine control commands are sent to the machine actuators to perform human-machine collaboration and adversarial control.

9. The human-machine collaboration and adversarial control optimization method for complex dynamic environments according to claim 6, characterized in that, The determination of the judgment threshold includes: Residual spectral features are extracted from multiple historical time windows corresponding to the cross-source residual correlation matrix before the event to form a historical singular value sequence. The pre-event residual strength benchmark is determined based on the preset quantile value of the historical singular value sequence; An adaptive judgment threshold is generated based on the pre-event residual strength benchmark and the preset safety margin; The singular values ​​corresponding to the current cross-source residual increment matrix are compared with the adaptive decision threshold, and the common residual directions that reach the adaptive decision threshold are determined as post-event confounding features.

10. The human-machine collaboration and adversarial control optimization method for complex dynamic environments according to claim 1, characterized in that, Also includes S8: From the point of the adversarial event anchor, the update of the human-machine closed-loop baseline path model using post-event data is stopped, and the parameter freeze state is entered. The path credibility of the current data is calculated based on the path attribution coefficient vector and the state reconstruction error. When the time window for a consecutive preset number of joint mutation indicators is lower than the recovery threshold and the attribution coefficient corresponding to post-event confounding factors is lower than the preset limit, the parameter freeze state is lifted. Select samples whose path credibility reaches the update threshold from the data after the parameter freeze is lifted, and use them as credible update samples. The baseline path parameter matrix is ​​updated incrementally using trusted update samples to suppress the contamination of parameters in the human-machine closed-loop baseline path model by post-event mixed data.

Images