A fermentation abnormality identification method based on temperature and humidity sequence

Abstract

The application discloses a kind of based on temperature and humidity sequence's fermentation anomaly identification method, specifically related to fermentation anomaly identification technical field, including obtaining temperature and humidity data in fermentation process, the temperature and humidity value of each time is recorded according to fermentation batch, monitoring position and time label, forms temperature and humidity original sequence;Eliminate missing value and error value in temperature and humidity original sequence, resample according to uniform time interval;According to adjacent time temperature and humidity value variation quantity, calculate variation intensity index, and compare with preset variation intensity index threshold, identify stable analysis section, form temperature and humidity reference sequence;By regarding temperature and humidity sequence as the system output under the action of microbial community hierarchical structure, and in stable analysis section, temperature and humidity sequence is executed driving force separation and driving force structure parameter identification, to drive force structure mutation as abnormal basis, realize accurate identification to early instability state of liquor fermentation.

Description

Technical Field

[0001] This invention relates to the field of fermentation anomaly identification technology, and more specifically, to a fermentation anomaly identification method based on temperature and humidity sequences. Background Technology

[0002] In the solid-state fermentation process of baijiu, the industry generally understands the temperature and humidity sequence as a continuous physical change curve that characterizes the release of fermentation heat and the state of water migration. The mainstream technology relies on the rising, plateauing, and falling trends of temperature and humidity values, as well as threshold deviations, rate of change, or morphological inflection points to determine whether the fermentation process is abnormal. However, the real fermentation system is not a homogeneous reaction. The fermentation dynamics are composed of multiple functional microbial communities, which form a hierarchical synergistic and energy transfer organizational structure within the grain layer. The dominant bacteria are responsible for driving the metabolic rhythm, while the subordinate microbial communities perform maintenance reactions. There is a stable structural control relationship between different levels. The temperature and humidity curves do not exist independently, but are macroscopic outputs formed synchronously driven by this internal organizational structure. When the hierarchical structure of the microbial community is stable, the temperature and humidity sequence shows a continuous evolution driven by a single central force. When the functional division of microorganisms is disrupted, the dominant bacteria decline, or the cooperative transmission pathway is interrupted, the driving source changes from a single dominant source to multiple weak independent sources, and the temperature and humidity changes become superimposed responses, losing consistency and structural control, leading to the internal system entering an early stage of instability. Since the sudden change in the driving source usually occurs when the single point value is still within the normal range, the temperature and humidity curve still presents a stable, reasonable appearance without abnormal signals, thus making it impossible for identification methods based on thresholds, morphology or prediction models to discover that the internal structure of the system has collapsed. It is evident that existing technologies mistakenly regard temperature and humidity sequences as direct representations of fermentation status, ignoring their true identity as system outputs driven by the microbial community hierarchy. This leads to a fundamental misdefinition of the abnormal identification targets, causing existing methods to completely fail in the most critical early instability stage. Summary of the Invention

[0003] To overcome the above-mentioned defects of the prior art, embodiments of the present invention provide a fermentation anomaly identification method based on temperature and humidity sequences. By treating the temperature and humidity sequences as the system output under the action of the microbial community hierarchy, and performing driving force separation and driving force structure parameter identification on the temperature and humidity sequences within the stable analysis segment, the method uses the abrupt change in the driving force structure as the basis for anomaly identification, thereby achieving accurate identification of the early unstable state of baijiu fermentation.

[0004] To achieve the above objectives, the present invention provides the following technical solution: a method for identifying fermentation anomalies based on temperature and humidity sequences, comprising: S1. Acquire temperature and humidity data during the fermentation process, and record the temperature and humidity values ​​at each moment according to the fermentation batch, monitoring location and time label to form the original temperature and humidity sequence; S2. Remove missing and erroneous values ​​from the original temperature and humidity sequence, and resample according to a uniform time interval; calculate the change intensity index based on the change in temperature and humidity values ​​at adjacent times, and compare it with the preset change intensity index threshold to identify stable analysis segments and form a temperature and humidity baseline sequence. S3. For the temperature and humidity reference sequence, extract the temperature and humidity subsequences within each window according to the preset time window, use the blind source separation algorithm to decompose the temperature and humidity subsequences within each time window to obtain the driving force subsequence, and splice the driving force subsequences corresponding to each time window to form a set of driving force sequences. S4. Calculate the energy magnitude, peak position, and duration based on the sampled values ​​of the driving force sequence within the current time window; identify the number of dominant and subordinate driving forces within the current time window by comparing the energy magnitude, peak position, and duration of each driving force sequence, and calculate the cooperative duration between each driving force based on the high-energy overlap duration to obtain the driving force structure parameter sequence. S5. Compare the sequence of driving force structure parameters corresponding to all time windows in each fermentation batch with the preset threshold for the number of dominant driving forces and the threshold for the duration of synergy, identify the driving force structure mutation window, and determine the corresponding fermentation batch as fermentation abnormal.

[0005] In a preferred embodiment, in S1, before the start of baijiu fermentation, a fermentation batch number is assigned to each fermentation batch, a collection location number is assigned to each monitoring location, and the correspondence between the collection location number and the corresponding monitoring location in the fermentation pit is recorded. During the fermentation process of baijiu, temperature and humidity data at each monitoring location are collected at preset sampling time intervals. The temperature and humidity data include the temperature and humidity values ​​at the corresponding time. The temperature and humidity values ​​at each time are bound to the corresponding fermentation batch number, monitoring location number and collection time label and written into the data record table. The records in the data recording table are sorted according to the fermentation batch number, the collection location number, and the time label to generate a temperature and humidity time series for each fermentation batch and each collection location. All fermentation batches and the corresponding temperature and humidity time series of all collection locations were collected according to the fermentation batch number to form the original temperature and humidity sequence.

[0006] In a preferred embodiment, in S2, the temperature value, humidity value and corresponding time tag in the original temperature and humidity sequence are organized into a temperature and humidity record sequence sorted by time tag. Each temperature and humidity record in the temperature and humidity record sequence includes the temperature value and humidity value corresponding to a single time tag. Traverse the temperature and humidity record sequence, check the integrity of the temperature and humidity values ​​for each temperature and humidity record. When a null value or an unresolved abnormal mark is detected, mark the corresponding temperature and humidity record as a missing record and delete it from the temperature and humidity record sequence to obtain the first temperature and humidity record sequence after removing the missing records. For each temperature and humidity record in the first temperature and humidity record sequence, a validity check is performed on the temperature and humidity values. When the temperature value is lower than the preset lower limit threshold or higher than the preset upper limit threshold, or the humidity value is lower than the preset lower limit threshold or higher than the preset upper limit threshold, the corresponding temperature and humidity record is marked as an erroneous record and deleted, thus obtaining the second temperature and humidity record sequence.

[0007] In a preferred embodiment, S2 further includes generating a resampling time axis within the time range of the second temperature and humidity recording sequence according to a preset resampling time interval. For time tags on the resampling time axis that already contain temperature and humidity records, the original temperature and humidity values ​​are retained. For time tags on the resampling time axis that lack temperature and humidity records, the two most recent temperature and humidity records before and after the corresponding time tag are searched. The total number of resampling time steps between the previous and next temperature and humidity recording times is calculated. The current time tag is determined as the current step number relative to the previous temperature and humidity record. The difference between the temperature value of the previous and next temperature and humidity records is multiplied by the ratio of the current step number to the total number of steps to obtain the temperature increment. The temperature increment is added to the temperature value of the previous temperature and humidity record to obtain the interpolated temperature value of the current time tag. The interpolated humidity value of the current time tag is calculated in the same way, thus forming a resampling temperature and humidity sequence with a uniform time interval. For each pair of adjacent time tags in the resampled temperature and humidity sequence, the sum of the absolute values ​​of the temperature difference and humidity difference between the current time tag and the previous time tag is used as the change intensity index corresponding to the current time tag. The change intensity index is compared with the preset change intensity index threshold. When the change intensity index corresponding to at least two consecutive time tags is less than the preset change intensity index threshold, the corresponding time interval is marked as a stable analysis segment, and the temperature and humidity records corresponding to all stable analysis segments are summarized into a temperature and humidity baseline sequence. Otherwise, the corresponding time interval is marked as an unstable analysis segment.

[0008] In a preferred embodiment, in S3, the temperature and humidity reference sequence is divided according to a preset time window. For each time window, the temperature values ​​belonging to the current time window are arranged in the order of the time labels to form the temperature subsequence of the current time window, and the humidity values ​​belonging to the current time window are arranged in the order of the time labels to form the humidity subsequence of the current time window. The temperature and humidity subsequences within the current time window are aligned one by one according to the time labels. The temperature and humidity values ​​corresponding to each time label are combined to form a sample vector containing two components, thus forming a two-dimensional input sequence for the current time window. Perform blind source separation operation on the two-dimensional input sequence of the current time window, including: Iterate through all sample vectors within the current time window and calculate the number of samples as the total number of samples. Sum the temperature values ​​in all sample vectors and divide by the total number of samples to obtain the average temperature. Sum the humidity values ​​in all sample vectors and divide by the total number of samples to obtain the average humidity. Subtract the average temperature value and the average humidity value from the temperature value of each sample vector in the two-dimensional input sequence to obtain a centered sequence used to eliminate the overall offset; For each sample vector in the centered sequence, multiply the temperature value by itself and sum them. Divide the sum by the total number of samples to obtain the mean square value of the temperature change. For each sample vector in the centered sequence, multiply the humidity value by itself and sum them. Divide the sum by the total number of samples to obtain the mean square value of the humidity change. For each sample vector in the centered sequence, multiply the temperature value and humidity value by themselves and sum them. Divide the sum by the total number of samples to obtain the average coupling amount of the temperature and humidity changes. The mean square value of temperature change, the mean square value of humidity change, and the average coupling amount of temperature and humidity change are arranged into a two-row, two-column table of change correlation coefficients according to the positional relationship between temperature change and humidity change. The mean square value of temperature change in the first row and first column of the correlation coefficient table is compared with the mean square value of humidity change in the second row and second column. When the mean square value of temperature change is greater than the mean square value of humidity change, temperature change is taken as the main component and humidity change is taken as the secondary component; when the mean square value of humidity change is greater than the mean square value of temperature change, humidity change is taken as the main component and temperature change is taken as the secondary component. When constructing the directional parameter vector, the first component of the directional parameter vector is defined as representing the direction of temperature change, and the second component is defined as representing the direction of humidity change. Based on the main change component and the secondary change component, two directional parameter vectors containing the first and second components are constructed. The directional parameter vector representing the main change trend takes a value of one only at the component position corresponding to the main change component and a value of zero at the other component position. The directional parameter vector representing the secondary change trend takes a value of one only at the component position corresponding to the secondary change component and a value of zero at the other component position, thus forming two independent and perpendicular directional parameter vectors. Use the two direction parameter vectors as two initial separation direction vectors; When calculating the candidate driving force subsequence for the separation direction vector, the temperature value of each sample vector in the centered sequence is multiplied by the first component of the current separation direction vector, the humidity value of the corresponding sample vector is multiplied by the second component of the current separation direction vector, and the two products are added together to obtain the candidate driving force value corresponding to the sample. All candidate driving force values ​​are arranged in order of time label to form the candidate driving force subsequence. When updating the separation direction vector, each candidate driving force value in the candidate driving force subsequence is multiplied by the cube of the corresponding candidate driving force value, and the update amount is obtained by summing all the products. The update amount is then added to the first and second components of the current separation direction vector respectively. Multiply each component of the updated separation direction vector by itself and sum them to obtain the sum of squares. Take the square root of the sum of squares to obtain the vector length. Divide each component of the current separation direction vector by the corresponding vector length to obtain the normalized separation direction vector. During orthogonalization, the first component of the two normalized separation direction vectors is multiplied by the second component and the sum is obtained to get the similarity value. The similarity value is multiplied by the first and second components of another normalized separation direction vector to form the projection vector. The corresponding projection vector is subtracted component by component from the current normalized separation direction vector to obtain the separation direction vectors that are perpendicular to each other. When the component difference between two consecutive updates is lower than the preset convergence threshold, the separation direction vector that has been normalized and orthogonalized is taken as the converged separation direction vector, and each converged separation direction vector is multiplied by the centering sequence to obtain multiple independent driving force subsequences, and an output channel number is assigned to each driving force subsequence.

[0009] In a preferred embodiment, S3 further includes arranging each driving force subsequence obtained within the current time window according to the order of the time tags, constructing a driving force record entry containing the driving force value, time tag and output channel number, and writing all driving force record entries into the driving force record set corresponding to the current time window. The driving force record sets corresponding to all time windows are read sequentially according to the time order of the time windows, and the driving force record entries are grouped according to the output channel number; for driving force record entries belonging to adjacent time windows under the same output channel number, the driving force record entries are connected according to the time label order to form a continuous driving force sequence. All driving force sequences are grouped together according to the output channel number as the index key to form a driving force sequence set, which contains multiple driving force sequences that can be traversed sequentially in time.

[0010] In a preferred embodiment, in S4, for each driving force sequence within the current time window, the driving force sampling values ​​are traversed in the order of the time labels, and all driving force sampling values ​​are multiplied by themselves and summed to obtain the cumulative driving force energy value. The cumulative driving force energy value is divided by the number of sampling points within the current time window to obtain the driving force energy magnitude. For each driving force sequence within the current time window, find the sampling point with the upper limit value in the driving force sampling value, record the time label corresponding to the sampling point as the main peak position, take the value of the corresponding sampling point as the main peak value, and traverse the adjacent sampling points forward and backward on both sides of the main peak value to find the continuous sampling point interval where the driving force sampling value is greater than or equal to half of the main peak value, and take the time length of the corresponding continuous sampling point interval as the duration. All driving force sequences within the current time window are sorted in descending order of driving force energy. The driving force sequence ranked first is marked as the dominant driving force sequence. When the energy difference between the energy of multiple driving force sequences and the energy of the dominant driving force sequence is less than a preset energy difference threshold, the corresponding multiple driving force sequences are marked as the dominant driving force sequence, and the remaining driving force sequences are marked as subordinate driving force sequences, thereby determining the number of dominant driving forces and subordinate driving forces.

[0011] In a preferred embodiment, S4 further includes, for each driving force sequence within the current time window, traversing each sampling point in the corresponding driving force sequence, when the driving force sampling value is greater than a preset energy threshold, marking the time tag of the corresponding sampling point as a high-energy time tag, and arranging all the high-energy time tags in chronological order to form a high-energy time tag sequence of the corresponding driving force sequence. For each pair of different driving force sequences within the current time window, read the high-energy time tag sequences corresponding to the two driving force sequences respectively. Take each time tag in the high-energy time tag sequence of the first driving force sequence as the current comparison time tag. In the high-energy time tag sequence of the second driving force sequence, check whether there is a time tag that is the same as the current comparison time tag. When the same time tag is detected, increment the overlap count value by one. After completing the traversal of all high-energy time tags of the first driving force sequence, the overlap count value is multiplied by the preset sampling time interval to obtain the high-energy overlap duration of the corresponding driving force sequence in the current time window; and with each driving force sequence as the benchmark, the high-energy overlap durations corresponding to other driving force sequences are accumulated to obtain the cooperative duration of the corresponding driving force sequence in the current time window. For each driving force sequence within the current time window, the output channel number, driving force energy magnitude, main peak position, duration, dominant driving force or subordinate driving force mark, and cooperative duration corresponding to the driving force sequence are sequentially written into the same driving force structure parameter record, and all driving force structure parameter records are stored in the driving force structure parameter record set corresponding to the current time window. In the set of driving force structure parameter records, the driving force structure parameter records are first sorted according to the output channel number. Within the same output channel number, they are sorted from smallest to largest according to the time label corresponding to the main peak position. The sorted driving force structure parameter records are then connected sequentially to obtain the driving force structure parameter sequence of the current time window. The driving force structural parameter sequences obtained from multiple time windows are arranged in chronological order of the time windows to form a set of driving force structural parameter sequences.

[0012] In a preferred embodiment, in S5, for each fermentation batch, the driving force structure parameter sequence formed by all time windows of the corresponding fermentation batch is read, and the dominant driving force quantity field in each driving force structure parameter record in each time window is compared with the preset dominant driving force quantity threshold. When the value of the dominant driving force quantity field is greater than the preset dominant driving force quantity threshold, the corresponding time window is marked as a dominant driving force abnormal window. For each fermentation batch, read the driving force structure parameter sequence formed by all time windows of the corresponding fermentation batch, and compare the collaborative duration field of each driving force structure parameter record in each time window with the preset collaborative duration threshold. When the value corresponding to the collaborative duration field is less than the preset collaborative duration threshold, the corresponding time window is marked as a collaborative duration abnormal window. When processing each fermentation batch, read the time window number list marked as the dominant driving force anomalous window and the time window number list marked as the cooperative persistence anomalous window in the corresponding fermentation batch, and initialize the driving force structure mutation window list to an empty list. Take out each time window number in the time window number list of the dominant driving force anomaly window as the current time window number to be detected. Traverse each time window number in the time window number list of the cooperative continuous anomaly window. When a time window number with the same number as the current time window number to be detected is detected, add the current time window number to the driving force structure mutation window list. After traversing all time window numbers in the time window number list of the dominant driving force anomaly window, mark all time windows contained in the driving force structure mutation window list as driving force structure mutation windows. When the number of driving force structure mutation windows is not less than one, the corresponding fermentation batch is determined to be an abnormal batch of baijiu fermentation, and the time position corresponding to the driving force structure mutation window is output as the basis for locating the fermentation abnormality; otherwise, the corresponding fermentation batch is determined to be a normal batch of baijiu fermentation.

[0013] The technical effects and advantages of this invention are as follows: This scheme treats temperature and humidity sequences as external projections of the microbial community hierarchy. By blind source separation, it solves the driving force subsequence and identifies changes in dominant and subordinate driving forces. This allows for the accurate exposure of early instability caused by the shift of driving sources from a single dominant source to multiple sources and coordinated disintegration, fundamentally solving the problem that traditional numerical and trend methods cannot identify early structural anomalies. By performing centering, constructing a correlation coefficient table of changes, updating iterative directions, and orthogonalizing within the stable analysis segment, the scheme can separate multiple independent driving force subsequences from temperature and humidity sequences containing only two observation channels, revealing the true driving mode of internal organization structure and improving the observability of microbial dynamic changes. Based on the magnitude of driving force energy, the location of the main peak, the duration of the main peak, and the duration of the synergistic effect, a sequence of driving force structural parameters is constructed, which makes the internal functional hierarchy, the strong and weak master-slave relationship, and the synergistic integrity of each time window quantifiable. This provides a structured and traceable indicator system for systemic instability, rather than relying on single-point temperature and humidity values. By using the intersection of two conditions—the number threshold of dominant driving forces and the duration threshold of collaboration—anomalies are only marked when multiple dominant driving forces appear and the collaborative structure collapses simultaneously. This ensures that the judgment logic is consistent with the actual microbial structural mutation mechanism, significantly reducing false alarms and improving anomaly location accuracy. By eliminating missing and erroneous records, resampling via linear interpolation, and screening stable analysis segments, the integrity, consistency, and comparability of the temperature and humidity baseline sequence are ensured from the source. This avoids external noise interference with subsequent driving force extraction and structure determination, making the entire process reliable and feasible in an industrial data acquisition environment. Attached Figure Description

[0014] Figure 1 This is a flowchart outlining the method steps of the present invention. Detailed Implementation

[0015] 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.

[0016] Refer to the instruction manual appendix Figure 1 An embodiment of the present invention provides a fermentation anomaly identification method based on temperature and humidity sequences, comprising: S1. Acquire temperature and humidity data collected in chronological order from multiple monitoring locations during the fermentation process of baijiu. Record the temperature and humidity values ​​at each moment according to the fermentation batch, monitoring location, and time label to form the original temperature and humidity sequence for driving force structure analysis. The temperature and humidity data include the temperature and humidity values ​​at the corresponding moment. S2. Remove missing and erroneous values ​​from the original temperature and humidity sequence, and resample the original temperature and humidity sequence at a uniform time interval. During the removal process, calculate the change intensity index based on the change in temperature and humidity values ​​at adjacent times, and compare it with the preset change intensity index threshold to identify stable analysis segments and form a temperature and humidity baseline sequence. S3. For the temperature and humidity reference sequence, extract the temperature subsequence and humidity subsequence within each window according to the preset time window. Use the blind source separation algorithm to decompose the temperature subsequence and humidity subsequence within each time window to obtain multiple independent driving force subsequences. Then, splice the driving force subsequences corresponding to each time window to form a set of driving force sequences. S4. Calculate the energy magnitude, peak position, and duration based on the sampled values ​​of the driving force sequence within the current time window; based on the comparison of energy magnitude, peak position, and duration between each driving force sequence, identify the number of dominant driving forces and subordinate driving forces within the current time window; and calculate the cooperative duration between each driving force based on the overlap length of the time intervals in different driving force sequences where the energy is greater than the preset energy threshold, thereby obtaining the driving force structure parameter sequence for the current time window. S5. Compare the sequence of driving force structure parameters corresponding to all time windows in each fermentation batch with the preset threshold for the number of dominant driving forces and the preset threshold for the duration of synergy, identify the driving force structure mutation window, and determine the corresponding fermentation batch as an abnormal fermentation of baijiu.

[0017] In S1, before the start of baijiu fermentation, a fermentation batch number is assigned to each fermentation batch, a collection location number is assigned to each monitoring location, and the correspondence between the collection location number and the corresponding monitoring location in the fermentation pit is recorded. During the fermentation of baijiu, temperature and humidity data at each monitoring location are collected at preset sampling time intervals. The temperature and humidity data include the temperature and humidity values ​​at the corresponding time. The temperature and humidity values ​​at each time are bound to the corresponding fermentation batch number, monitoring location number and collection time label and written into the data record table. This is used to uniformly store and manage the raw temperature and humidity data of baijiu fermentation during the collection phase. The records in the data recording table are sorted according to the fermentation batch number, collection location number, and time tag to generate a temperature and humidity time series for each fermentation batch and each collection location, which reflects the temperature and humidity changes at that location throughout the entire baijiu fermentation cycle. All fermentation batches and temperature and humidity time series corresponding to all collection locations were collected according to the fermentation batch number to form the original temperature and humidity sequence for driving force structure analysis. It should be noted that, in this scheme, a fermentation batch refers to a group of grain mash units that are fed into the same cellar or fermentation unit within the same opening time and fermented in solid-state baijiu according to the same process conditions. The fermentation batch number assigned to each fermentation batch is used to uniquely identify the baijiu fermentation process at the data level, ensuring that the temperature and humidity data collected from different monitoring locations for the same batch can be uniformly collected in subsequent analysis. The monitoring location refers to the specific physical location where temperature and humidity sensors are pre-selected and installed inside or on the surface of the baijiu fermentation pit, including but not limited to the upper, middle, and lower layers of the pit, as well as locations near the side walls of the pit. The acquisition location number assigned to each monitoring location is used to identify the monitoring location without directly recording the physical coordinates, and to record the correspondence between the acquisition location number and the corresponding monitoring location in the pit, so that the temperature and humidity distribution can be restored according to the acquisition location number in subsequent analysis. Temperature and humidity data refer to the temperature and humidity values ​​collected by temperature and humidity sensors at preset sampling time intervals during the fermentation process of baijiu. The temperature value is used to characterize the temperature state of the grain mash and the surrounding environment at the monitoring location, and the humidity value is used to characterize the relative humidity state of the air or grain layer at the monitoring location. Together, they constitute a quantitative description of the microenvironment of baijiu fermentation at the monitoring location. Time stamps refer to the time information recorded by the acquisition device or host computer system when collecting temperature and humidity data. They are used to sort multiple temperature and humidity data from the same fermentation batch and the same acquisition location by time and form a time series. In the specific implementation process, before the start of baijiu fermentation, each fermentation batch is assigned a fermentation batch number, and each monitoring location is assigned a collection location number, with the correspondence between the two recorded. During the baijiu fermentation process, temperature and humidity values ​​at each monitoring location are collected at preset sampling time intervals. Each collected temperature and humidity value is bound to the corresponding fermentation batch number, collection location number, and time label and written into a data record table. When the baijiu fermentation ends or the predetermined data processing cycle arrives, the records in the data record table are sorted according to the order of fermentation batch number, collection location number, and time label, generating a temperature and humidity time series for each fermentation batch and each collection location. The temperature and humidity time series corresponding to all collection locations under the same fermentation batch are aggregated to form the original temperature and humidity sequence for driving force structure analysis. This establishes the correlation between fermentation batch, monitoring location, and time series at the data structure level, providing a complete and consistent data foundation for subsequent driving force decomposition and baijiu fermentation anomaly identification based on temperature and humidity sequences.

[0018] In S2, the temperature values, humidity values, and corresponding time tags in the original temperature and humidity sequence are organized into a temperature and humidity record sequence sorted by time tags. Each temperature and humidity record in the temperature and humidity record sequence includes the temperature value and humidity value corresponding to a single time tag. The temperature and humidity record sequence is traversed, and the integrity of the temperature and humidity values ​​is checked for each record. When a null temperature or humidity value or an unresolved abnormal marker is detected, the corresponding temperature and humidity record is marked as a missing record and deleted from the temperature and humidity record sequence, resulting in the first temperature and humidity record sequence after removing the missing records. The null temperature or humidity value or unresolved abnormal marker refers to the invalid placeholder value returned by the acquisition device when the acquisition fails, the sensor is offline, or communication packets are lost. This includes empty strings, blank fields, out-of-limit padding values, or invalid data markers identified by a fixed encoding format, which are used to indicate that the temperature and humidity values ​​are unavailable at that sampling time. For each temperature and humidity record in the first temperature and humidity record sequence, a validity check is performed on the temperature and humidity values. When the temperature value is lower than the preset lower temperature threshold or higher than the preset upper temperature threshold, or the humidity value is lower than the preset lower humidity threshold or higher than the preset upper humidity threshold, the corresponding temperature and humidity record is marked as an erroneous record and deleted. This results in a second temperature and humidity record sequence containing only valid temperature and humidity values. This ensures that subsequent resampling interpolation processing is performed only on valid data, avoiding interference from erroneous records in the calculation of change intensity indicators and the identification of stable analysis segments.

[0019] S2 also includes generating a resampling time axis within the time range of the second temperature and humidity recording sequence according to a preset resampling time interval. For time tags on the resampling time axis that already contain temperature and humidity records, the original temperature and humidity values ​​are retained. For time tags on the resampling time axis that lack temperature and humidity records, the two most recent temperature and humidity records before and after the corresponding time tag are searched. The total number of resampling time steps between the previous and next temperature and humidity record times is calculated. The current time tag is determined as the current resampling step relative to the previous temperature and humidity record. The difference between the temperature value of the previous and next temperature and humidity records is multiplied by the ratio of the current step to the total number of steps to obtain the temperature increment. The temperature increment is added to the temperature value of the previous temperature and humidity record to obtain the interpolated temperature value of the current time tag. The interpolated humidity value of the current time tag is calculated in the same way, thus forming a resampling temperature and humidity sequence with a uniform time interval. The two most recent temperature and humidity records refer to the temperature and humidity record with a time tag earlier than the current time tag and closest to the current time tag in the second temperature and humidity record sequence, and the temperature and humidity record with a time tag later than the current time tag and closest to the current time tag, which are used to define the reference points before and after the interpolation calculation. For each pair of adjacent time tags in the resampled temperature and humidity sequence, the sum of the absolute values ​​of the temperature difference and humidity difference between the current time tag and the previous time tag is used as the change intensity index corresponding to the current time tag. The change intensity index is compared with the preset change intensity index threshold. When the change intensity index corresponding to at least two consecutive time tags is less than the preset change intensity index threshold, the corresponding time interval is marked as a stable analysis segment. The temperature and humidity records corresponding to all stable analysis segments are summarized into a temperature and humidity reference sequence, which is used as the temperature and humidity reference sequence input to the driving force structure analysis step to provide a reference for the temperature and humidity changes of baijiu fermentation under stable background conditions. Otherwise, the corresponding time interval is marked as an unstable analysis segment.

[0020] In S3, the temperature and humidity reference sequence is divided according to a preset time window. For each time window, the temperature values ​​belonging to the current time window are arranged in the order of the time labels to form the temperature subsequence of the current time window, and the humidity values ​​belonging to the current time window are arranged in the order of the time labels to form the humidity subsequence of the current time window. The temperature and humidity subsequences within the current time window are aligned one by one according to the time labels. The temperature and humidity values ​​corresponding to each time label are combined to form a sample vector containing two components, thus forming a two-dimensional input sequence for the current time window. Perform blind source separation operation on the two-dimensional input sequence of the current time window, including: Iterate through all sample vectors within the current time window and calculate the number of samples as the total number of samples. Sum the temperature values ​​in all sample vectors and divide by the total number of samples to obtain the average temperature. Sum the humidity values ​​in all sample vectors and divide by the total number of samples to obtain the average humidity. Subtract the average temperature value and the average humidity value from the temperature value of each sample vector in the two-dimensional input sequence to obtain a centered sequence used to eliminate the overall offset; For each sample vector in the centered sequence, multiply the temperature value by itself and sum them. Divide the sum by the total number of samples to obtain the mean square value of the temperature change. For each sample vector in the centered sequence, multiply the humidity value by itself and sum them. Divide the sum by the total number of samples to obtain the mean square value of the humidity change. For each sample vector in the centered sequence, multiply the temperature value and humidity value by themselves and sum them. Divide the sum by the total number of samples to obtain the average coupling amount of the temperature and humidity changes. The average square value of temperature change, the average square value of humidity change, and the average coupling amount of temperature and humidity change are arranged into a two-row, two-column correlation coefficient table according to the positional relationship between temperature change and humidity change. The two-row, two-column correlation coefficient table is a matrix structure formed by arranging the average square value of temperature change, the average square value of humidity change, and the average coupling amount of temperature and humidity change in the fermentation of liquor according to the positional relationship between temperature component and humidity component. The first row and first column are the average square value of temperature change, the second row and second column are the average square value of humidity change, and the first row and second column are the average coupling amount of temperature and humidity change, which is used to measure the linear coupling relationship between temperature and humidity change. The mean square value of temperature change in the first row and first column of the correlation coefficient table is compared with the mean square value of humidity change in the second row and second column. When the mean square value of temperature change is greater than the mean square value of humidity change, temperature change is taken as the main component and humidity change is taken as the secondary component; when the mean square value of humidity change is greater than the mean square value of temperature change, humidity change is taken as the main component and temperature change is taken as the secondary component. When constructing the directional parameter vector for driving force separation calculation, the first component of the directional parameter vector is defined as representing the direction of temperature change, and the second component is defined as representing the direction of humidity change. Based on the main change component and the secondary change component, two directional parameter vectors containing the first and second components are constructed. The directional parameter vector representing the main change trend takes a value of one only at the component position corresponding to the main change component and a value of zero at the other component position. The directional parameter vector representing the secondary change trend takes a value of one only at the component position corresponding to the secondary change component and a value of zero at the other component position, thus forming two independent and mutually perpendicular directional parameter vectors. The two directional parameter vectors are used as two initial separation directional vectors, and these vectors are used as the first separation directional vectors for the calculation and iterative update of the candidate driving force subsequence, providing the initial directional basis for subsequent calculation of driving force values ​​through sample projection and execution of separation directional update operations. When calculating the candidate driving force subsequence for the separation direction vector, the temperature value of each sample vector in the centered sequence is multiplied by the first component of the current separation direction vector, the humidity value of the corresponding sample vector is multiplied by the second component of the current separation direction vector, and the two products are added together to obtain the candidate driving force value corresponding to the sample. All candidate driving force values ​​are arranged in order of time label to form the candidate driving force subsequence. When updating the separation direction vector, each candidate driving force value in the candidate driving force subsequence is multiplied by the cube of the corresponding candidate driving force value, and the update amount is obtained by summing all the products. The update amount is then added to the first and second components of the current separation direction vector respectively. Multiply each component of the updated separation direction vector by itself and sum them to obtain the sum of squares. Take the square root of the sum of squares to obtain the vector length. Divide each component of the current separation direction vector by the corresponding vector length to obtain the normalized separation direction vector. During orthogonalization, the first component of the two normalized separation direction vectors is multiplied by the second component and the sum is obtained to get the similarity value. The similarity value is multiplied by the first and second components of another normalized separation direction vector to form the projection vector. The corresponding projection vector is subtracted component by component from the current normalized separation direction vector to obtain the separation direction vectors that are perpendicular to each other. When the component difference between two consecutive updates is lower than the preset convergence threshold, the separation direction vector that has been normalized and orthogonalized is taken as the converged separation direction vector, and each converged separation direction vector is multiplied by the centering sequence to obtain multiple independent driving force subsequences, and an output channel number is assigned to each driving force subsequence.

[0021] S3 also includes arranging each driving force subsequence obtained in the current time window according to the order of the time tags, constructing a driving force record entry containing the driving force value, time tag and output channel number, and writing all driving force record entries into the driving force record set corresponding to the current time window, which is used to store the driving force record entries under all output channel numbers in the current time window; The driving force record sets corresponding to all time windows are read sequentially according to the time order of the time windows, and the driving force record entries are grouped according to the output channel number; for driving force record entries belonging to adjacent time windows under the same output channel number, the driving force record entries are connected according to the time label order to form a driving force sequence arranged continuously throughout the entire analysis time range; All driving force sequences are grouped together according to the output channel number as the index key to form a driving force sequence set. The driving force sequence set contains multiple driving force sequences that can be traversed sequentially in time. The driving force sequence set is used as input data for subsequent driving force structure parameter identification calculations and baijiu fermentation anomaly determination calculations.

[0022] In S4, for each driving force sequence within the current time window, the driving force sample values ​​are traversed in the order of the time labels. All driving force sample values ​​are multiplied by themselves and summed to obtain the cumulative driving force energy value. The cumulative driving force energy value is divided by the number of sampling points within the current time window to obtain the driving force energy magnitude. Here, the driving force sample value refers to the driving force value corresponding to each time label of the driving force sequence within the current time window. For each driving force sequence within the current time window, find the sampling point with the upper limit value in the driving force sampling value, record the time label corresponding to the sampling point as the main peak position, take the value of the corresponding sampling point as the main peak value, and traverse the adjacent sampling points forward and backward on both sides of the main peak value to find the continuous sampling point interval where the driving force sampling value is greater than or equal to half of the main peak value, and take the time length of the corresponding continuous sampling point interval as the duration. All driving force sequences within the current time window are sorted in descending order of driving force energy. The driving force sequence ranked first is marked as the dominant driving force sequence. When the energy difference between the energy of multiple driving force sequences and the energy of the dominant driving force sequence is less than a preset energy difference threshold, the corresponding multiple driving force sequences are marked as dominant driving force sequences, and the remaining driving force sequences are marked as subordinate driving force sequences. This determines the number of dominant and subordinate driving forces, which is used to distinguish the master and slave roles of different driving forces in the driving force structure parameters within the current time window, providing a basis for calculating the duration of collaboration and constructing the driving force structure parameter record.

[0023] S4 also includes, for each driving force sequence within the current time window, traversing each sampling point in the corresponding driving force sequence, when the driving force sampling value is greater than the preset energy threshold, marking the time label of the corresponding sampling point as a high-energy time label, and arranging all the high-energy time labels in chronological order to form a high-energy time label sequence of the corresponding driving force sequence; For each pair of different driving force sequences within the current time window, the high-energy time tag sequences corresponding to the two driving force sequences are read. Each time tag in the high-energy time tag sequence of the first driving force sequence is used as the current comparison time tag. The high-energy time tag sequence of the second driving force sequence is searched for a time tag that is the same as the current comparison time tag. When a time tag that is the same is detected, the overlap count is incremented by one. It should be noted that the overlap count is initialized to zero before processing each pair of driving force sequences. After traversing all high-energy time tags of the first driving force sequence, the overlap count is multiplied by the preset sampling time interval to obtain the high-energy overlap duration of the corresponding driving force sequence in the current time window. Then, based on each driving force sequence, the high-energy overlap durations corresponding to other driving force sequences are accumulated to obtain the cooperative duration of the corresponding driving force sequence in the current time window. This is used to characterize the cumulative time that the driving force and other driving forces maintain a high-energy state simultaneously in the high-energy interval, providing a numerical basis for the cooperative duration field in the driving force structure parameter record. For each driving force sequence within the current time window, the output channel number, driving force energy magnitude, main peak position, duration, dominant driving force or subordinate driving force mark, and cooperative duration corresponding to the driving force sequence are sequentially written into the same driving force structure parameter record, and all driving force structure parameter records are stored in the driving force structure parameter record set corresponding to the current time window. In the set of driving force structure parameter records, the driving force structure parameter records are first sorted according to the output channel number. Within the same output channel number, they are sorted from smallest to largest according to the time label corresponding to the main peak position. The sorted driving force structure parameter records are then connected sequentially to obtain the driving force structure parameter sequence of the current time window. The driving force structure parameter sequences obtained from multiple time windows are arranged in chronological order of the time windows to form a set of driving force structure parameter sequences for identifying driving force structure mutations and determining abnormalities in baijiu fermentation.

[0024] In S5, for each fermentation batch, the driving force structure parameter sequence formed by all time windows of the corresponding fermentation batch is read. The dominant driving force quantity field in each driving force structure parameter record in each time window is compared with the preset dominant driving force quantity threshold. When the value of the dominant driving force quantity field is greater than the preset dominant driving force quantity threshold, the corresponding time window is marked as a dominant driving force abnormal window. The dominant driving force quantity field refers to the number of dominant driving force sequences that are statistically obtained and written into the driving force structure parameter record during the construction of the driving force structure parameter record for the current time window. It is used to characterize the number of entries marked as dominant driving force sequences within the time window. For each fermentation batch, the driving force structure parameter sequence formed by all time windows of the corresponding fermentation batch is read. The cooperative duration field in each driving force structure parameter record in each time window is compared with the preset cooperative duration threshold. When the value corresponding to the cooperative duration field is less than the preset cooperative duration threshold, the corresponding time window is marked as a cooperative duration abnormal window. The cooperative duration field refers to the value obtained and written into the driving force structure parameter record during the calculation of the driving force cooperative duration. It is used to characterize the cumulative length of time that the current driving force sequence is in a high-energy state at the same time as other driving force sequences within the current time window. When processing each fermentation batch, the time window number lists marked as dominant driving force anomaly windows and cooperative persistence anomaly windows in the corresponding fermentation batch are read respectively. The driving force structure mutation window list is initialized to an empty list. Initializing the list to an empty list means creating a new data list to store time window numbers for the driving force structure mutation window. No time window numbers are written to this data list when it is created, so that the number of elements in the data list is zero. This allows the corresponding time window numbers to be appended sequentially when a time window that meets the conditions is detected later. The time window number list refers to the set of numbers assigned to each time window in the process of constructing the driving force structure parameter sequence to uniquely identify the time window according to its chronological order within the Baijiu fermentation cycle. The time window numbers are sequentially incremented in chronological order. Each time window number has a one-to-one correspondence with the time range of the driving force structure parameter records within the corresponding time window. It is used as an index reference when performing anomaly window judgment and window location across time windows. Take out each time window number in the time window number list of the dominant driving force anomaly window as the current time window number to be detected. Traverse each time window number in the time window number list of the cooperative continuous anomaly window. When a time window number with the same number as the current time window number to be detected is detected, add the current time window number to the driving force structure mutation window list. After traversing all time window numbers in the time window number list of the dominant driving force anomaly window, mark all time windows contained in the driving force structure mutation window list as driving force structure mutation windows. When the number of driving force structure mutation windows is not less than one, the corresponding fermentation batch is determined to be an abnormal batch of baijiu fermentation, and the time position corresponding to the driving force structure mutation window is output as the basis for locating the fermentation abnormality; otherwise, the corresponding fermentation batch is determined to be a normal batch of baijiu fermentation.

[0025] The method for determining the values ​​of each preset threshold needs to be explained as follows: The preset lower temperature threshold and the preset upper temperature threshold are determined based on the allowed fermentation temperature range in the target liquor fermentation process document. Alternatively, the lower and upper quantile values ​​of the normal fermentation batch temperature records under the same yeast, the same cellar type, and the same season can be read as the value boundaries. The preset lower humidity threshold and the preset upper humidity threshold are determined according to the humidity range allowed in the target liquor fermentation process document. Alternatively, the lower and upper quantile values ​​of the humidity record of a normal fermentation batch can be read as the value boundaries. The preset threshold for the intensity of change index is determined based on the distribution of the intensity of change index of adjacent time tags in the temperature and humidity sequence of the normal fermentation batch resampling. Any value between the median value and the upper quartile value of the intensity of change index of the normal fermentation batch is taken as the screening boundary of the stable analysis segment. The preset convergence threshold is determined based on the numerical precision of the separation direction vector. It is less than one percent of the smaller of the standard deviations of the temperature and humidity components in the centered sequence. It is used to determine whether two consecutive separation direction updates should be stopped. The preset energy difference threshold is determined based on the distribution of energy differences between different driving force sequences within a normal fermentation batch. The lower quartile value of the energy difference between the first driving force and the other driving forces within the same time window of a normal fermentation batch is taken as the boundary for determining parallel dominance. The preset energy threshold is determined based on the distribution of the squared values ​​of the driving force sampling values ​​within the current fermentation batch. The average energy value of the driving force sampling points within the same time window is taken as the high-energy time label screening boundary, or the average energy value of the driving force sampling points within the corresponding time window of the normal fermentation batch is taken as the screening boundary. The preset threshold for the number of dominant driving forces is determined based on the statistical results of the number of dominant driving forces in each time window of the normal fermentation batch. The maximum occurrence value of the number of dominant driving forces in the normal fermentation batch is taken as the threshold. When the number of dominant driving forces in the batch to be identified is greater than the threshold, it indicates that the driving force structure has changed from a single dominant state to a multi-dominant state. The preset threshold for the duration of collaboration is determined based on the statistical results of the duration of collaboration of each driving force sequence in the normal fermentation batch. The lower quartile value of the duration of collaboration of the normal fermentation batch is taken as the threshold. When the duration of collaboration of the batch to be identified is less than the threshold, it indicates that the high-energy overlap relationship between the driving forces is weakened. Therefore, it can be seen that the aforementioned preset thresholds can be determined based on the target liquor fermentation process documents, the statistical distribution of historical normal fermentation batches, the numerical accuracy of the data acquisition equipment, and conventional process experience in the field of liquor fermentation. Those skilled in the art can obtain the corresponding threshold values ​​according to the target cellar type, yeast type, seasonal conditions, and sampling interval.

[0026] It is important to clarify that, including but not limited to: traditional methods treat temperature and humidity curves as direct indicators of fermentation status, focusing only on numerical values, trends, and threshold deviations, while ignoring the fact that the temperature and humidity sequence is actually a macroscopic output driven by the internal microbial community hierarchy; when the dominant microbial community declines and synergistic relationships are disrupted, the driving source changes from a single dominant source to multiple weakly independent sources, and temperature and humidity changes become superimposed responses that lose consistency. During this stage, the values ​​are often still within the normal range, causing existing threshold or trend models to completely fail in the most critical early instability stage; based on this understanding, the process of forming the solution is: no longer treating temperature and humidity as the result itself, but as a projection of the driving force structure, focusing on how to extract the changes in the driving force structure from the temperature and humidity sequence to form the overall link from S1 to S5; First, at the data organization level of S1, the solution uses three dimensions—fermentation batch number, collection location number, and time tag—to organize the originally mixed temperature and humidity collection records in the baijiu fermentation pit into an original temperature and humidity sequence that is ordered by time, located by space, and distinguishable by batch. The fermentation batch is defined as a group of grain mash units that are fermented in solid-state baijiu under the same opening time, the same pit or fermentation unit, and under the same process conditions. The collection location corresponds to the specific location of the upper, middle, lower layers and side walls of the pit. The time tag is recorded by the collection device or host computer for sorting. Secondly, before identifying anomalies in S2, the reliability of the temperature and humidity data itself must be addressed. This solution does not directly cram the collected values ​​into the algorithm. Instead, it first forms a temperature and humidity record sequence according to time tags, and then checks each temperature and humidity value to see if it is empty, or if it is an abnormal placeholder value due to acquisition failure or communication packet loss. Missing records are deleted first. Then, the remaining records are validated, and obviously out-of-bounds erroneous records are removed to obtain a second temperature and humidity record sequence containing only valid temperature and humidity values. The reason for doing this is that sensor disconnection, poor cable contact, and upper computer communication interruption are very common in real-world cellar environments. If the data is not cleaned first, the subsequent change intensity indicators and driving force decomposition will be severely interfered with by these garbage points. Next, instead of simply using the original sampling time interval, this solution introduces a resampling time axis. Within the time range of the second sequence, a time axis is generated according to a uniform resampling interval, and missing time points are filled in using linear interpolation: For each time tag without a record, the two most recent temperature and humidity records before and after it are found, and the number of resampling time steps between these two records is counted. Then, the position of the current time tag in these resampling steps is determined, and the ratio of this current step number to the total number of steps is used to split the temperature and humidity differences between the two records before and after, thereby obtaining the interpolated temperature and humidity values ​​for the current time point. This resampling logic only uses addition, subtraction, comparison, and proportional allocation, and can be implemented by any industrial host computer. Finally, in order to perform driving force decomposition only in the relatively stable range of microbial community structure, this scheme calculates the sum of the absolute values ​​of temperature difference and humidity difference for each pair of adjacent time tags in the resampled temperature and humidity sequence as a change intensity index. When the change intensity of at least two consecutive time tags is lower than the preset threshold, the corresponding time interval is marked as a stable analysis segment, and the temperature and humidity records corresponding to all stable analysis segments are summarized into a temperature and humidity baseline sequence. Building upon this, S3 addresses how to separate the hidden driving force structure from the temperature and humidity baseline sequence. First, the temperature and humidity baseline sequence is divided into preset time windows. Within each window, temperature and humidity values ​​are arranged in chronological order, forming temperature and humidity subsequences. Then, at each time label, the temperature and humidity values ​​are concatenated into a sample vector containing two components. Thus, the sequence of sample vectors within a time window constitutes a two-dimensional input sequence. This embodies a core assumption: the temperature and humidity at each moment are two observation channels acting together by the same set of internal driving forces. Only by placing them in a two-dimensional space can true source separation be achieved. Next, the blind source separation operation is broken down into a series of fully feasible operational procedures: First, iterate through all sample vectors within the window to calculate the average temperature and average humidity. Then, for each sample, subtract the corresponding average from the temperature and humidity to obtain a centered sequence. This step eliminates the overall bias, ensuring that what is seen later is the changing structure rather than the absolute level. Then, based on the centered sequence, multiply and sum the temperature values ​​of each sample, divide by the total number of samples to obtain the average square value of the temperature change. Similarly, obtain the average square value of the humidity change. Then, multiply and sum the temperature and humidity values ​​of each sample, divide by the total number of samples to obtain the average coupling amount of the temperature and humidity changes. Use these three quantities to construct a two-row, two-column table of change correlation coefficients according to the positional relationship of the temperature and humidity components. This table is essentially a centralized expression of the autocorrelation of temperature changes, the autocorrelation of humidity changes, and the linear coupling between temperature and humidity. Then, by comparing the magnitudes of the mean square of temperature change and the mean square of humidity change in the table, the dimension with the more drastic change is taken as the principal component and the other dimension as the secondary component. Based on this, a direction parameter vector is constructed: the first component of the direction parameter vector corresponds specifically to the direction of temperature change, and the second component corresponds to the direction of humidity change. The principal component takes a value of one at the position corresponding to the principal change and zero at the other position. The secondary component takes a value of one at the position corresponding to the secondary change and zero at the other position. In this way, two independent and mutually perpendicular initial separation direction vectors are obtained. In each subsequent iteration, these direction vectors and the centered sequence are used to perform very specific calculations: For each sample, the temperature value is multiplied by the first component of the current direction vector, the humidity value is multiplied by the second component, and then the results are added together to obtain the candidate driving force value for that sample. All candidate driving force values ​​are arranged in chronological order to form a candidate driving force subsequence. Then, the sum of the products of each candidate driving force value and its own cube value is used as the update amount, which is added to each component of the current direction vector. The updated vector is then normalized to ensure that the magnitude of the direction vector is one. For two normalized separation direction vectors, the similarity is obtained by multiplying and summing the components. This similarity is then used to construct a projection vector, and this projection is subtracted from the current vector component by component to achieve orthogonality between the two direction vectors. When the change of each component between two consecutive updates is less than the preset convergence threshold, the separation direction is considered to have converged. The same projection calculation is then performed on the centered sequence using the converged direction vector to obtain multiple independent driving force subsequences, and an output channel number is assigned to each subsequence. Between time windows, this scheme further uses driving force record entries and driving force record sets to splice the driving force subsequences scattered in different time windows under the same output channel number into a continuous driving force sequence along the entire analysis time range. All driving force sequences are then grouped into a driving force sequence set according to the output channel number, which serves as the input for subsequent structural parameter identification and baijiu fermentation anomaly determination.

[0027] This highlights that the instability of the microbial community hierarchy is often not an instantaneous change within a window, but rather a slow accumulation over time. Therefore, it is essential to ensure that the temporal trajectory of the same underlying driving force is continuously traceable between windows. After obtaining the set of driving force sequences, S4 first iterates through the sampled values ​​of each driving force sequence in chronological order within each time window, multiplies each driving force sampled value by itself and accumulates them, and then divides them by the number of sampling points in the window to define the energy magnitude; by finding the maximum value in the sampled values ​​and its corresponding time label as the main peak position, and by finding a continuous interval on the left and right of the main peak where the driving force sampled value is greater than or equal to half of the main peak to define the duration, the three together characterize the intensity and duration of each driving force within the window; Then, the driving force sequence is arranged from largest to smallest energy. The one with the highest energy is marked as the dominant driving force sequence. When the energy difference between multiple driving forces and this one is less than a preset energy difference threshold, they are all marked as dominant driving forces, and the rest are marked as subordinate driving forces. In this way, the question of whether the system is dominated by one force or torn apart by multiple forces of the same level is transformed into a change in the number of dominant driving forces. Meanwhile, to characterize the synergistic relationship between different driving forces, this scheme selects time tags with energy values ​​greater than a preset energy threshold at sampling points for each driving force sequence, forming a high-energy time tag sequence. For each pair of different driving force sequences, the high-energy time tags of the first sequence are traversed, and the same time tags are searched in the high-energy time tag sequence of the second sequence. Each time the same tag is found, the overlap count is incremented by one. Finally, the overlap count is multiplied by the sampling time interval to obtain the high-energy overlap duration of the two driving forces within the window. Then, taking each driving force sequence as a benchmark, the high-energy overlap duration between it and all other driving forces is accumulated to obtain the synergistic duration of the driving force within the current window. This synergistic duration field corresponds to the question mentioned in the background technology about whether the hierarchical synergy of the microbial community and energy transfer are still intact: if there are many dominant driving forces but the high-energy overlap between them is very short, it indicates that the system has already slid from a single dominant, strong synergy to a multi-source, weak synergy unstable state. Finally, the fields such as output channel number, driving force energy magnitude, main peak position, duration, dominant or subordinate marker, and synergistic duration are written into the driving force structural parameter record set. Within the set, the parameters are first sorted by output channel number and then by main peak time to form the driving force structural parameter sequence for the current time window. The structural parameter sequences of each time window are then arranged in chronological order to form the driving force structural parameter sequence set, providing a two-dimensional perspective of time and structure for the final batch-level judgment. S5 abstracts this set of structural parameters into a simple but very powerful judgment logic: For each fermentation batch, it reads the sequence of driving force structural parameters for all time windows of that batch, compares the dominant driving force quantity field with the preset dominant driving force quantity threshold window by window, and if the dominant driving force quantity exceeds the threshold in a certain window, the window is recorded as a dominant driving force abnormal window; similarly, it compares the synergistic duration field with the preset synergistic duration threshold window by window, and if the synergistic duration is lower than the threshold in a certain window, the window is recorded as a synergistic duration abnormal window; Then, within each batch, the intersection of these two types of abnormal window lists is calculated: the time window numbers are sequentially retrieved from the dominant driving force abnormal window list, and the same number is checked in the cooperative persistence abnormal window list. If the same number exists, it is added to the driving force structure mutation window list. After traversal, all time windows contained in the mutation window list are marked as driving force structure mutation windows. If the number of mutation windows is not less than one, the fermentation batch is judged as an abnormal batch of baijiu fermentation, and the time positions corresponding to these windows are output as the basis for anomaly location. Otherwise, it is judged as a normal batch of baijiu fermentation. The core here is to take the two conditions of the driving source changing from a single dominant source to multiple sources and the cooperative structure disintegration as the premise for anomaly judgment. Only when the two conditions are met simultaneously in the same time window is it considered that the internal driving force structure has undergone a substantial mutation. This corresponds to the analysis of the instability mechanism in the background, and the threshold and list intersection method ensures clarity, controllability and adjustability in implementation. In summary, this approach does not modify existing threshold judgments from the outset. Instead, it redesigns the entire chain based on the view that temperature and humidity sequences are external projections of the microbial community hierarchy. It establishes a unified spatial, temporal, and batch-specific data framework using fermentation batch numbers, sampling location numbers, and time tags. It constructs a temperature and humidity baseline sequence through anomaly removal and resampling. It obtains the driving force subsequence by performing a clearly defined blind source separation step within the stable analysis period and constructs a set of driving force sequences covering the entire fermentation cycle. It then uses structural parameters such as energy, peak position, duration, number of dominant driving forces, and synergistic duration to characterize the internal organizational state. Finally, it identifies the driving force structural mutation window through a combination of dominant driving force number thresholds and synergistic duration thresholds, enabling the localizable determination of early instability states in baijiu fermentation batches.

[0028] The above description is merely a preferred embodiment of the present invention and is not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the protection scope of the present invention.

Claims

1. A method for identifying fermentation anomalies based on temperature and humidity sequences, characterized in that, include: S1. Acquire temperature and humidity data during the fermentation process, and record the temperature and humidity values ​​at each moment according to the fermentation batch, monitoring location and time label to form the original temperature and humidity sequence; S2. Remove missing and erroneous values ​​from the original temperature and humidity sequence, and resample according to a uniform time interval; calculate the change intensity index based on the change in temperature and humidity values ​​at adjacent times, and compare it with the preset change intensity index threshold to identify stable analysis segments and form a temperature and humidity baseline sequence. S3. For the temperature and humidity reference sequence, extract the temperature and humidity subsequences within each window according to the preset time window, use the blind source separation algorithm to decompose the temperature and humidity subsequences within each time window to obtain the driving force subsequence, and splice the driving force subsequences corresponding to each time window to form a set of driving force sequences. Specifically, the temperature and humidity two-dimensional sequence is centered, a correlation coefficient table is constructed to determine the main and secondary change components, and an initial separation direction vector is constructed. The temperature and humidity values ​​of the centered sequence are multiplied by the components of the current separation direction vector, and the two products are added to obtain the candidate driving force value corresponding to the sample. Each candidate driving force value in the candidate driving force subsequence is multiplied by the cube of the corresponding candidate driving force value, and the sum of all products is obtained to obtain the update value. The updated separation direction vector is normalized and orthogonalized, and iterated until the component difference is lower than the threshold and convergence occurs. The converged direction vector is multiplied by the centered sequence to output multiple independent driving force subsequences and assign channel numbers. S4. Calculate the energy magnitude, main peak position, and duration based on the sampled values ​​of the driving force sequence within the current time window; Based on the comparison of energy magnitude, peak position and duration among the driving force sequences, the number of dominant and subordinate driving forces within the current time window is identified, and the cooperative duration among the driving forces is calculated based on the high-energy overlap duration, thus obtaining the driving force structure parameter sequence. S5. Compare the sequence of driving force structure parameters corresponding to all time windows in each fermentation batch with the preset threshold for the number of dominant driving forces and the threshold for the duration of synergy, identify the driving force structure mutation window, and determine the corresponding fermentation batch as fermentation abnormal.

2. The fermentation anomaly identification method based on temperature and humidity sequence according to claim 1, characterized in that: In S1, before the start of baijiu fermentation, a fermentation batch number is assigned to each fermentation batch, a collection location number is assigned to each monitoring location, and the correspondence between the collection location number and the corresponding monitoring location in the fermentation pit is recorded. During the fermentation process of baijiu, temperature and humidity data at each monitoring location are collected at preset sampling time intervals. The temperature and humidity data include the temperature and humidity values ​​at the corresponding time. The temperature and humidity values ​​at each time are bound to the corresponding fermentation batch number, monitoring location number and collection time label and written into the data record table. The records in the data recording table are sorted according to the fermentation batch number, the collection location number, and the time label to generate a temperature and humidity time series for each fermentation batch and each collection location. All fermentation batches and the corresponding temperature and humidity time series of all collection locations were collected according to the fermentation batch number to form the original temperature and humidity sequence.

3. The fermentation anomaly identification method based on temperature and humidity sequence according to claim 2, characterized in that: In S2, the temperature values, humidity values, and corresponding time tags in the original temperature and humidity sequence are organized into a temperature and humidity record sequence sorted by time tags. Each temperature and humidity record in the temperature and humidity record sequence includes the temperature value and humidity value corresponding to a single time tag. Traverse the temperature and humidity record sequence, check the integrity of the temperature and humidity values ​​for each temperature and humidity record. When a null value or an unresolved abnormal mark is detected, mark the corresponding temperature and humidity record as a missing record and delete it from the temperature and humidity record sequence to obtain the first temperature and humidity record sequence after removing the missing records. For each temperature and humidity record in the first temperature and humidity record sequence, a validity check is performed on the temperature and humidity values. When the temperature value is lower than the preset lower limit threshold or higher than the preset upper limit threshold, or the humidity value is lower than the preset lower limit threshold or higher than the preset upper limit threshold, the corresponding temperature and humidity record is marked as an erroneous record and deleted, thus obtaining the second temperature and humidity record sequence.

4. The fermentation anomaly identification method based on temperature and humidity sequence according to claim 3, characterized in that: S2 also includes generating a resampling time axis within the time range of the second temperature and humidity recording sequence according to a preset resampling time interval. For time tags on the resampling time axis that already contain temperature and humidity records, the original temperature and humidity values ​​are retained. For time tags on the resampling time axis that lack temperature and humidity records, the two most recent temperature and humidity records before and after the corresponding time tag are searched. The total number of resampling time steps between the previous and next temperature and humidity recording times is calculated. The current time tag is determined as the current step number relative to the previous temperature and humidity record. The difference between the temperature value of the previous and next temperature and humidity records is multiplied by the ratio of the current step number to the total number of steps to obtain the temperature increment. The temperature increment is added to the temperature value of the previous temperature and humidity record to obtain the interpolated temperature value of the current time tag. The interpolated humidity value for the current time tag is calculated using the same method, thus forming a resampled temperature and humidity sequence with a uniform time interval; For each pair of adjacent time tags in the resampled temperature and humidity sequence, the sum of the absolute values ​​of the temperature difference and humidity difference between the current time tag and the previous time tag is used as the change intensity index corresponding to the current time tag. The change intensity index is compared with the preset change intensity index threshold. When the change intensity index corresponding to at least two consecutive time tags is less than the preset change intensity index threshold, the corresponding time interval is marked as a stable analysis segment, and the temperature and humidity records corresponding to all stable analysis segments are summarized into a temperature and humidity baseline sequence. Otherwise, the corresponding time interval is marked as an unstable analysis segment.

5. The fermentation anomaly identification method based on temperature and humidity sequence according to claim 4, characterized in that: In S3, the temperature and humidity reference sequence is divided according to a preset time window. For each time window, the temperature values ​​belonging to the current time window are arranged in the order of the time labels to form the temperature subsequence of the current time window, and the humidity values ​​belonging to the current time window are arranged in the order of the time labels to form the humidity subsequence of the current time window. The temperature and humidity subsequences within the current time window are aligned one by one according to the time labels. The temperature and humidity values ​​corresponding to each time label are combined to form a sample vector containing two components, thus forming a two-dimensional input sequence for the current time window. Perform blind source separation operation on the two-dimensional input sequence of the current time window, including: Iterate through all sample vectors within the current time window and calculate the number of samples as the total number of samples. Sum the temperature values ​​in all sample vectors and divide by the total number of samples to obtain the average temperature. Sum the humidity values ​​in all sample vectors and divide by the total number of samples to obtain the average humidity. Subtract the average temperature value and the average humidity value from the temperature value of each sample vector in the two-dimensional input sequence to obtain a centered sequence used to eliminate the overall offset; For each sample vector in the centered sequence, multiply the temperature value by itself and sum them. Divide the sum by the total number of samples to obtain the mean square value of the temperature change. For each sample vector in the centered sequence, multiply the humidity value by itself and sum them. Divide the sum by the total number of samples to obtain the mean square value of the humidity change. For each sample vector in the centered sequence, multiply the temperature value and humidity value by themselves and sum them. Divide the sum by the total number of samples to obtain the average coupling amount of the temperature and humidity changes. The mean square value of temperature change, the mean square value of humidity change, and the average coupling amount of temperature and humidity change are arranged into a two-row, two-column table of change correlation coefficients according to the positional relationship between temperature change and humidity change. The mean square value of temperature change in the first row and first column of the correlation coefficient table is compared with the mean square value of humidity change in the second row and second column. When the mean square value of temperature change is greater than the mean square value of humidity change, temperature change is taken as the main component and humidity change is taken as the secondary component; when the mean square value of humidity change is greater than the mean square value of temperature change, humidity change is taken as the main component and temperature change is taken as the secondary component. When constructing the directional parameter vector, the first component of the directional parameter vector is defined as representing the direction of temperature change, and the second component is defined as representing the direction of humidity change. Based on the main change component and the secondary change component, two directional parameter vectors containing the first and second components are constructed. The directional parameter vector representing the main change trend takes a value of one only at the component position corresponding to the main change component and a value of zero at the other component position. The directional parameter vector representing the secondary change trend takes a value of one only at the component position corresponding to the secondary change component and a value of zero at the other component position, thus forming two independent and perpendicular directional parameter vectors. Use the two direction parameter vectors as two initial separation direction vectors; When calculating the candidate driving force subsequence for the separation direction vector, the temperature value of each sample vector in the centered sequence is multiplied by the first component of the current separation direction vector, the humidity value of the corresponding sample vector is multiplied by the second component of the current separation direction vector, and the two products are added together to obtain the candidate driving force value corresponding to the sample. All candidate driving force values ​​are arranged in order of time label to form the candidate driving force subsequence. When updating the separation direction vector, each candidate driving force value in the candidate driving force subsequence is multiplied by the cube of the corresponding candidate driving force value, and the update amount is obtained by summing all the products. The update amount is then added to the first and second components of the current separation direction vector respectively. Multiply each component of the updated separation direction vector by itself and sum them to obtain the sum of squares. Take the square root of the sum of squares to obtain the vector length. Divide each component of the current separation direction vector by the corresponding vector length to obtain the normalized separation direction vector. During orthogonalization, the first component of the two normalized separation direction vectors is multiplied by the second component and the sum is obtained to get the similarity value. The similarity value is multiplied by the first and second components of another normalized separation direction vector to form the projection vector. The corresponding projection vector is subtracted component by component from the current normalized separation direction vector to obtain the separation direction vectors that are perpendicular to each other. When the component difference between two consecutive updates is lower than the preset convergence threshold, the separation direction vector that has been normalized and orthogonalized is taken as the converged separation direction vector, and each converged separation direction vector is multiplied by the centering sequence to obtain multiple independent driving force subsequences, and an output channel number is assigned to each driving force subsequence.

6. The fermentation anomaly identification method based on temperature and humidity sequence according to claim 5, characterized in that: S3 also includes arranging each driving force subsequence obtained within the current time window according to the order of the time tags, constructing a driving force record entry containing the driving force value, time tag and output channel number, and writing all driving force record entries into the driving force record set corresponding to the current time window; The driving force record sets corresponding to all time windows are read sequentially according to the time order of the time windows, and the driving force record entries are grouped according to the output channel number; for driving force record entries belonging to adjacent time windows under the same output channel number, the driving force record entries are connected according to the time label order to form a continuous driving force sequence. All driving force sequences are grouped together according to the output channel number as the index key to form a driving force sequence set, which contains multiple driving force sequences that can be traversed sequentially in time.

7. The fermentation anomaly identification method based on temperature and humidity sequence according to claim 6, characterized in that: In S4, for each driving force sequence within the current time window, the driving force sampling values ​​are traversed in the order of the time label. All driving force sampling values ​​are multiplied by themselves and summed to obtain the cumulative driving force energy value. The cumulative driving force energy value is divided by the number of sampling points within the current time window to obtain the driving force energy magnitude. For each driving force sequence within the current time window, find the sampling point with the upper limit value in the driving force sampling value, record the time label corresponding to the sampling point as the main peak position, take the value of the corresponding sampling point as the main peak value, and traverse the adjacent sampling points forward and backward on both sides of the main peak value to find the continuous sampling point interval where the driving force sampling value is greater than or equal to half of the main peak value, and take the time length of the corresponding continuous sampling point interval as the duration. All driving force sequences within the current time window are sorted in descending order of driving force energy. The driving force sequence ranked first is marked as the dominant driving force sequence. When the energy difference between the energy of multiple driving force sequences and the energy of the dominant driving force sequence is less than a preset energy difference threshold, the corresponding multiple driving force sequences are marked as the dominant driving force sequence, and the remaining driving force sequences are marked as subordinate driving force sequences, thereby determining the number of dominant driving forces and subordinate driving forces.

8. The fermentation anomaly identification method based on temperature and humidity sequence according to claim 7, characterized in that: S4 also includes, for each driving force sequence within the current time window, traversing each sampling point in the corresponding driving force sequence, when the driving force sampling value is greater than the preset energy threshold, marking the time label of the corresponding sampling point as a high-energy time label, and arranging all the high-energy time labels in chronological order to form a high-energy time label sequence of the corresponding driving force sequence; For each pair of different driving force sequences within the current time window, read the high-energy time tag sequences corresponding to the two driving force sequences respectively. Take each time tag in the high-energy time tag sequence of the first driving force sequence as the current comparison time tag. In the high-energy time tag sequence of the second driving force sequence, check whether there is a time tag that is the same as the current comparison time tag. When the same time tag is detected, increment the overlap count value by one. After completing the traversal of all high-energy time tags of the first driving force sequence, the overlap count value is multiplied by the preset sampling time interval to obtain the high-energy overlap duration of the corresponding driving force sequence in the current time window; and with each driving force sequence as the benchmark, the high-energy overlap durations corresponding to other driving force sequences are accumulated to obtain the cooperative duration of the corresponding driving force sequence in the current time window. For each driving force sequence within the current time window, the output channel number, driving force energy magnitude, main peak position, duration, dominant driving force or subordinate driving force mark, and cooperative duration corresponding to the driving force sequence are sequentially written into the same driving force structure parameter record, and all driving force structure parameter records are stored in the driving force structure parameter record set corresponding to the current time window. In the set of driving force structure parameter records, the driving force structure parameter records are first sorted according to the output channel number. Within the same output channel number, they are sorted from smallest to largest according to the time label corresponding to the main peak position. The sorted driving force structure parameter records are then connected sequentially to obtain the driving force structure parameter sequence of the current time window. The driving force structural parameter sequences obtained from multiple time windows are arranged in chronological order of the time windows to form a set of driving force structural parameter sequences.

9. The fermentation anomaly identification method based on temperature and humidity sequence according to claim 8, characterized in that: In S5, for each fermentation batch, the driving force structure parameter sequence formed by all time windows of the corresponding fermentation batch is read. The dominant driving force quantity field in each driving force structure parameter record in each time window is compared with the preset dominant driving force quantity threshold. When the value of the dominant driving force quantity field is greater than the preset dominant driving force quantity threshold, the corresponding time window is marked as a dominant driving force abnormal window. For each fermentation batch, read the driving force structure parameter sequence formed by all time windows of the corresponding fermentation batch, and compare the collaborative duration field of each driving force structure parameter record in each time window with the preset collaborative duration threshold. When the value corresponding to the collaborative duration field is less than the preset collaborative duration threshold, the corresponding time window is marked as a collaborative duration abnormal window. When processing each fermentation batch, read the time window number list marked as the dominant driving force anomalous window and the time window number list marked as the cooperative persistence anomalous window in the corresponding fermentation batch, and initialize the driving force structure mutation window list to an empty list. Take out each time window number in the time window number list of the dominant driving force anomaly window as the current time window number to be detected. Traverse each time window number in the time window number list of the cooperative continuous anomaly window. When a time window number with the same number as the current time window number to be detected is detected, add the current time window number to the driving force structure mutation window list. After traversing all time window numbers in the time window number list of the dominant driving force anomaly window, mark all time windows contained in the driving force structure mutation window list as driving force structure mutation windows. When the number of driving force structure mutation windows is not less than one, the corresponding fermentation batch is determined to be an abnormal batch of baijiu fermentation, and the time position corresponding to the driving force structure mutation window is output as the basis for locating the fermentation abnormality; otherwise, the corresponding fermentation batch is determined to be a normal batch of baijiu fermentation.

Images