A cloud-based approach to genome-wide genetic diversity analysis based on workflow automation

Abstract

This invention discloses a cloud-based method for whole-genome genetic diversity analysis based on automated workflow, specifically in the field of gene diversity analysis. It addresses the problem of pruning imbalance by dividing read length segments based on sequence beginning and ending features and constructing quantile tags. Quality control and pruning steps achieve accurate characterization of the original read length structure. A mapping matrix is ​​used to correlate quality distribution with segment features, and a two-parameter risk mapping model using variable point potential energy and trough symmetry measures quantifies the truncation sensitivity coefficient, enabling accurate identification of hidden over-pruned segments. A sensitive index set drives the update of a slow-release pruning strategy, restoring stability in alignment coverage and variant integrity after replacing static thresholds. Combining coverage curves and variant distribution data, parameter caching is dynamically updated, forming a self-correcting process. A collaborative data-driven workflow is constructed to effectively smooth out abnormal peaks in the diversity analysis process and optimize cloud-based parallel stability and resource utilization efficiency.

Description

Technical Field

[0001] This invention relates to the field of genetic diversity analysis, and more specifically, to a cloud-based method for whole-genome genetic diversity analysis based on workflow automation. Background Technology

[0002] In cross-laboratory scenarios, a large number of raw sequencing files are sequentially processed by FastQC and Trimmomatic via a cloud workflow to complete quality checks and trimming. The front-end interface directly incorporates early short-read sequencing experience thresholds into the workflow, but these thresholds are not linked to the real-time read length distribution, thus lacking the ability to adaptively adjust to differences in sample type and read length.

[0003] However, when uploaded samples contain excessively long read segments, the static threshold continuously triggers the pruning window, indiscriminately truncating functional sections. This results in coverage gaps during the comparison phase, missing fields in the mutation detection output, and ultimately, abnormal peaks on the population diversity curve, leading to incorrect classifications in the population structure analysis. Furthermore, the threshold is not dynamically adjusted with real-time read length statistics, and the logs are fragmented by parallel threads, making it difficult to trace the cause and impacting both process stability and result reliability.

[0004] To address the aforementioned problems, a technical solution is provided. Summary of the Invention

[0005] To overcome the aforementioned deficiencies of existing technologies, embodiments of the present invention provide a cloud-based method for whole-genome genetic diversity analysis based on automated workflow. This method divides read length segments based on sequence beginning and ending features and constructs quantile tags. Quality control and pruning steps achieve accurate characterization of the original read length structure. A mapping matrix is ​​used to correlate quality distribution with segment features. A two-parameter risk mapping model using variable point potential energy and trough symmetry measures quantifies the truncation sensitivity coefficient, enabling accurate identification of hidden over-truncated segments. A sensitive index set drives the update of a slow-release pruning strategy, restoring stability in alignment coverage and variant integrity after replacing static thresholds. Combined with coverage curves and variant distribution data, parameter caching is dynamically updated, forming a self-correcting process. A collaborative data-driven workflow is constructed to effectively smooth out abnormal peaks in the diversity analysis process, optimize cloud-based parallel stability and resource utilization efficiency, thereby solving the problems mentioned in the background technology.

[0006] To achieve the above objectives, the present invention provides the following technical solution:

[0007] S1: By analyzing the read length distribution and the first and last base quality characteristics of the original sequencing files, the read lengths are divided into different segments and corresponding quantile tags are generated;

[0008] S2: By calling the quality histogram distribution output by FastQC and combining it with quantile labels, a high- and low-quality read segment mapping matrix is ​​constructed to mark the quality stability of each read length segment;

[0009] S3: Determine the trigger frequency of sliding window pruning in different read length intervals based on the mapping matrix, calculate the segment trigger dispersion and coverage gap balance, synthesize the truncation sensitivity coefficient through the preset pruning risk mapping model, and extract the sensitive area index set;

[0010] S4: Backtrack the sensitive region index set to the original read length segment position, generate a pruning and release adjustment table and inject the pruning task configuration to replace the static threshold;

[0011] S5: Perform the adjusted pruning and compare the output sequence, extract the coverage curve and variation distribution, verify the coverage stability of sensitive segments and the integrity of the variation field, and then write the results to the parameter cache table.

[0012] In a preferred embodiment, step S1 includes the following:

[0013] Raw sequencing files in FASTQ format were collected, and read length distribution was analyzed to determine minimum length, maximum length, and distribution characteristics. Reads were segmented into intervals based on the length and quality characteristics of the start and end bases, and quantile labels were assigned to each interval according to length quantiles and quality categories.

[0014] In a preferred embodiment, step S2 includes the following:

[0015] The raw sequencing files were analyzed using FastQC to generate a quality histogram report, quality characteristic indicators were extracted, a high- and low-quality read mapping matrix was constructed using quantile tags, and stability scores were calculated by setting thresholds and nonlinear fusion to mark the quality stability of each read segment. Finally, the mapping matrix and stability tags were stored as structured data.

[0016] In a preferred embodiment, step S2 further includes the following:

[0017] Quality characteristics include average quality score, quality score fluctuation range, GC content deviation, and repetitive sequence ratio. The average quality score is obtained by calculating the average quality score of all bases within the read segment. The quality score fluctuation range is obtained by calculating the difference between the maximum and minimum quality scores within the read segment. The GC content deviation is determined by calculating the percentage absolute deviation between the GC content of the read segment and the global average GC content. The repetitive sequence ratio is obtained by calculating the proportion of repetitive sequences in the read segment to the total sequence.

[0018] In a preferred embodiment, step S3 includes the following:

[0019] By analyzing the mapping matrix, the trigger frequency of sliding window pruning in different read length intervals is determined, the segment trigger dispersion and coverage gap balance are calculated, and the truncation sensitivity coefficient is synthesized using a preset pruning risk mapping model. Read length segments with truncation sensitivity coefficients higher than a set threshold are extracted to form a sensitive region index set.

[0020] In a preferred embodiment, step S3 further includes the following:

[0021] The processing logic for segment trigger dispersion is as follows: the distribution of trigger events is obtained by using the sliding window trigger frequency sequence of each read length segment in the statistical mapping matrix. Then, the frequency sequence is subjected to variable point potential energy analysis to identify key points of frequency change. The dispersion index is obtained by calculating the fluctuation variance, thereby quantitatively describing the uneven distribution of clipping triggers.

[0022] In a preferred embodiment, step S3 further includes the following:

[0023] The processing logic for coverage gap balance is as follows: extract the segment coverage surface from the comparison curves before and after trimming, analyze the change characteristics of the coverage range, and then perform a peak-valley symmetry measurement on the segment coverage surface. By calculating the degree of symmetry between the peaks and valleys, the gap balance index is obtained, thereby assessing the impact of the trimming operation on the coverage integrity.

[0024] In a preferred embodiment, step S4 includes the following:

[0025] A sensitive region location table is generated by tracing the sensitive region index set back to the start and end positions of sensitive reads in the original sequencing file. A release factor is calculated using the quality stability score and truncation sensitivity coefficient, and a dynamic threshold is generated by combining this with the original static threshold to construct a pruning release adjustment table. The dynamic threshold in the pruning release adjustment table is integrated with the original static threshold into the pruning task configuration, and the pruning tool configuration file is updated to apply the dynamic threshold to sensitive reads and the original static threshold to non-sensitive reads.

[0026] In a preferred embodiment, step S5 includes the following:

[0027] The original sequencing file is trimmed using the dynamic and static threshold adjustment parameters configured in the trimming task, generating a trimmed sequencing file. The trimmed sequencing file is then aligned with the reference genome using the alignment tool BWA, generating an alignment file. Coverage depth sequences of sensitive regions are extracted from the alignment file, and the coverage curve fluctuation coefficient is calculated to quantify the degree of fluctuation of the coverage curve. The coverage curve fluctuation coefficient is compared with a preset stability threshold to determine the coverage stability of sensitive regions.

[0028] In a preferred embodiment, step S5 further includes the following:

[0029] The GATK mutation detection tool is used to compare and detect mutations in the file, generate mutation detection results, and check the integrity of the mutation fields in sensitive sections. Based on the integrity check results of the mutation fields, the integrity of the mutation fields in sensitive sections is determined. The coverage stability and mutation integrity status of sensitive sections are recorded in the parameter cache table, and the status is marked according to the verification results.

[0030] The technical effects and advantages of this invention's cloud-based method for analyzing whole-genome genetic diversity based on workflow automation are as follows:

[0031] By segmenting read length regions based on sequence beginning and end features and constructing quantile labels, the quality control and pruning steps achieve accurate characterization of the original read length structure. A mapping matrix is ​​used to correlate quality distribution with segment features, and a two-parameter risk mapping model measuring variable point potential energy and trough symmetry is used to quantify the truncation sensitivity coefficient, enabling accurate identification of hidden over-truncated segments. A sensitive index set drives the gradual pruning strategy update, restoring stability in alignment coverage and mutation integrity after replacing the static threshold. Combined with coverage curves and mutation distribution data, parameter caching is dynamically updated, forming a self-correcting process. From initial read length feature capture, quality mapping, risk quantification to pruning backtracking verification, a collaborative data-driven workflow is constructed, effectively smoothing out abnormal peaks in the diversity analysis process, reducing the need for manual intervention, and optimizing cloud parallel stability and resource utilization efficiency. Attached Figure Description

[0032] Figure 1 This is a flowchart illustrating the cloud-based method for whole-genome genetic diversity analysis based on workflow automation, as described in this invention. Detailed Implementation

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

[0034] Example 1: Figure 1 This invention presents a cloud-based method for whole-genome genetic diversity analysis based on workflow automation, including:

[0035] S1: By analyzing the read length distribution and the first and last base quality characteristics of the original sequencing files, the read lengths are divided into different segments and corresponding quantile tags are generated.

[0036] S2: By calling the quality histogram distribution output by FastQC and combining it with quantile labels, a high- and low-quality read segment mapping matrix is ​​constructed to mark the quality stability of each read length segment.

[0037] S3: Determine the trigger frequency of sliding window pruning in different read length intervals based on the mapping matrix, calculate the segment trigger dispersion and coverage gap balance, synthesize the truncation sensitivity coefficient through the preset pruning risk mapping model, and extract the sensitive area index set.

[0038] S4: Backtrack the sensitive region index set to the original read length segment position, generate a pruning and release adjustment table and inject the pruning task configuration to replace the static threshold.

[0039] S5: Perform the adjusted pruning and compare the output sequence, extract the coverage curve and variation distribution, verify the coverage stability of sensitive segments and the integrity of the variation field, and then write the results to the parameter cache table.

[0040] In whole-genome genetic diversity analysis, quality control and preprocessing of sequencing data are crucial for ensuring the accuracy of subsequent alignment and variant detection. In cross-laboratory scenarios, a large number of raw sequencing files are processed through cloud workflows. Traditional methods rely on fixed quality checks and trimming thresholds, such as using FastQC to assess quality and Trimmomatic to perform trimming. However, such static thresholds cannot adapt to the diversity of sample types and read length distributions. Especially when dealing with ultra-long reads, fixed thresholds may lead to over-trimming, resulting in truncated functional regions, coverage gaps during alignment, and thus abnormal results in variant detection and population diversity analysis. To address this issue, this invention proposes a cloud-based whole-genome genetic diversity analysis method based on workflow automation. By dynamically adjusting the trimming threshold to adapt to read length differences, the stability and reliability of the analysis are improved. Step S1, as the starting point of the entire process, is responsible for preliminary analysis of the raw sequencing data, dividing read length segments based on sequence beginning and ending characteristics and generating quantile tags, providing accurate basic data for subsequent quality assessment and trimming.

[0041] The purpose of step S1 is to divide the read length segments according to the continuous base combination pattern at the beginning and end of the sequence after collecting and uploading the original sequencing file of the sample, and to generate quantile tags simultaneously as a benchmark for subsequent quality evaluation.

[0042] Step S1 includes the following:

[0043] S1.1, Collect raw sequencing files:

[0044] Raw sequencing data files are obtained from user-uploaded samples, typically in FASTQ format. A FASTQ file contains the nucleotide sequence for each read and its corresponding quality fraction. The nucleotide sequence consists of bases, and the quality fraction consists of a quality value for each base. The quality value reflects the reliability of identifying that base during sequencing and is usually expressed numerically.

[0045] S1.2, Statistical Read Length Distribution:

[0046] The lengths of all reads in the original sequencing files are statistically analyzed to calculate the length range of the read set, including the minimum length value, the maximum length value, and the overall characteristics of the length distribution. The overall characteristics of the length distribution refer to the central tendency and dispersion of the read lengths. The central tendency is represented by the mean or median, and the dispersion is represented by the range of length values ​​or the standard deviation.

[0047] The purpose of statistically analyzing read length distribution is to comprehensively understand the variations in read length within a sample. By determining the minimum and maximum length values ​​and distribution characteristics, the overall properties of the read length set can be revealed, which helps identify typical length ranges and outliers, thereby improving the rationality and relevance of segmentation.

[0048] S1.3, Divide the reading length into segments:

[0049] Based on the read length and the quality characteristics of consecutive start and end bases, the read set is divided into multiple read segments. The specific process includes the following two steps:

[0050] First, define the quality characteristics of consecutive bases at the beginning and end of the read. Select a fixed number of bases at the start and end positions of the read, for example, 10 bases each. Calculate the average mass fraction of the first 10 bases, called the head mass fraction, and the average mass fraction of the last 10 bases, called the tail mass fraction.

[0051] Then, the partitioning method is implemented. The first step is to divide the reads into preliminary intervals based on length. The length dividing points are determined according to the read length distribution, dividing the read set into short read length segments, medium read length segments, and long read length segments. The dividing points are determined based on the quantiles of the length distribution, such as tertiles. The second step is to calculate the difference between the tail quality score and the head quality score for each read within each length interval; this is called the head-tail quality score difference.

[0052] Based on the magnitude of the difference between the first and last quality scores, the segments are further subdivided: if the difference is close to zero, it is classified as a high-quality sub-segment; if the difference deviates from zero but is within a preset range, it is classified as a medium-quality sub-segment; if the difference exceeds the preset range, it is classified as a low-quality sub-segment. This ultimately forms multiple read length segments, such as "short read length - high-quality segment" or "medium read length - low-quality segment".

[0053] Segmenting read lengths by combining length and initial / final quality features enables detailed classification of read lengths, accurately characterizing the structural differences within a set of reads. This method considers not only length variations but also quality distribution, identifying read length groups with different characteristics, thus providing targeted support for subsequent quality assessment and data processing. Compared to single-dimensional segmentation, multi-dimensional feature classification significantly improves the accuracy of the analysis.

[0054] S1.4, Generate quantile labels:

[0055] Quantile labels are generated for each read segment as a benchmark for subsequent quality assessment. The specific method is as follows: First, all reads in the original sequencing file are sorted by length from smallest to largest, and length quantile intervals are calculated. For example, the read length is divided into four quantile intervals, corresponding to the length ranges of 0% to 25%, 25% to 50%, 50% to 75%, and 75% to 100%, respectively. Then, the length range of each read segment is mapped to the corresponding quantile interval, and combined with its first and last quality characteristics, such as high quality, medium quality, or low quality, quantile labels are generated. Quantile labels are represented in the form of "length quantile interval - quality category," such as "first quantile interval - high quality" or "second quantile interval - low quality." The label design must ensure compatibility with the quality distribution output by subsequent quality analysis tools.

[0056] Unique labels are generated by linking read length segments to the global length distribution using a quantile method. The process first calculates the quantile intervals for all read lengths, then maps the length range and quality characteristics of each segment to the corresponding interval, forming structured labels. The calculation emphasizes the combination of global and local data to ensure that the labels reflect the relative position and quality characteristics of the read length segments.

[0057] Step S1 systematically processes sequencing data by collecting raw sequencing files, statistically analyzing read length distribution, dividing read length segments, and generating quantile tags. The read sets are classified into multiple segments based on length and beginning / end quality characteristics, and a standardized tag is assigned to each segment. This method comprehensively characterizes read length diversity, improves the accuracy and relevance of data analysis, and provides structured data support for subsequent processing.

[0058] Step S1 has already collected the raw sequencing files and divided the reads into segments based on the continuous base combination pattern at the beginning and end of the sequence, generating quantile tags to provide a benchmark for subsequent quality assessment. However, relying solely on the division of read segments and quantile tags cannot fully reveal the quality stability of each segment, especially when dealing with ultra-long reads, where static thresholds may lead to over-pruning, affecting the reliability of the analysis. Therefore, Step S2 needs to further utilize the quality distribution information output by FastQC, combined with quantile tags, to construct a high- and low-quality read mapping matrix and mark the quality stability of each read segment, providing accurate data support for the pruning frequency determination in Step S3.

[0059] Step S2 aims to call the quality histogram distributions output by FastQC, combine them with the quantile labels generated in step S1, construct a high- and low-quality read segment mapping matrix, and mark the quality stability of each read length segment in the matrix.

[0060] Step S2 includes the following:

[0061] S2.1, Call FastQC to output the quality histogram distribution:

[0062] FastQC analysis is performed on the acquired raw sequencing files to generate a quality histogram report. The quality histogram report includes the following information: base quality score distribution, recording statistical data on the quality score of each base position; sequence quality score distribution, providing the overall quality score for each read; GC content distribution, showing the GC content characteristics of the read segments; and repetitive sequence content, statistically analyzing the proportion of potential sequencing redundancy or contamination. After analysis, a report file containing the above statistical data is generated.

[0063] S2.2, Extract quality characteristic indicators:

[0064] The following quality characteristics are extracted from the quality histogram report to assess the quality status of each read segment: average quality score, obtained by calculating the average quality score of all bases within the read segment; quality score fluctuation range, obtained by calculating the difference between the maximum and minimum quality scores within the read segment; GC content deviation, determined by calculating the absolute percentage deviation between the GC content of the read segment and the global average GC content; and repetitive sequence ratio, obtained by calculating the proportion of repetitive sequences within the read segment to the total sequence.

[0065] S2.3, Construct a mapping matrix by combining quantile labels:

[0066] Using quantile labels, the quality characteristic indicators of each read length segment are associated with a matrix structure to construct a high-quality and low-quality read segment mapping matrix. The rows of the matrix correspond to quantile labels, representing different read length segments; the columns correspond to quality characteristic indicators, including average quality score, quality score fluctuation range, GC content deviation, and repetitive sequence ratio. Each element of the matrix is ​​the value of the corresponding read length segment on a specific quality characteristic indicator.

[0067] S2.4, Marking quality stability:

[0068] Based on the high-quality and low-quality read segment mapping matrix, the quality stability of each read length segment is quantitatively labeled, and the result is represented by a stability score. The specific steps are as follows: First, set thresholds for quality characteristic indicators. For example, the average quality score must be higher than 25, the quality score fluctuation range less than 12, the GC content deviation less than 6%, and the repetitive sequence ratio less than 2%. Then, determine whether each quality characteristic indicator meets the corresponding threshold; if it does, record it as 1, otherwise record it as 0. Next, calculate the stability score. For the average quality score that meets the conditions, take the natural logarithm. The quality score fluctuation range is processed using an exponential decay method (i.e., the value is decreased exponentially and the reciprocal is taken). The GC content deviation and the repetitive sequence ratio are processed using linear normalization (i.e., the value is mapped to the range of 0 to 1). Finally, add the processed values ​​and divide by 4 to obtain the stability score. Based on the stability score, read length segments are labeled into three categories: scores greater than or equal to 0.75 are high quality, scores between 0.45 and less than 0.75 are medium quality, and scores less than 0.45 are low quality.

[0069] S2.5, Output mapping matrix and labels:

[0070] The high- and low-quality read segment mapping matrix and the quality stability marker for each read length segment are stored as structured data for subsequent steps. The output data adopts a standardized format to ensure that it contains matrix elements and stability markers, facilitating subsequent steps to determine the pruning trigger frequency based on this data.

[0071] Step S2 generates a quality histogram distribution report by calling FastQC, extracts quality feature indicators, constructs a high- and low-quality read segment mapping matrix by combining quantile labels, and quantifies the quality stability score of each read length segment using a nonlinear fusion method, ultimately outputting structured data. This process achieves accurate quality stability labeling for the diverse characteristics of ultra-long read segment samples.

[0072] Step S2 calls FastQC to output the quality distribution and combines it with quantile labels to construct a high- and low-quality read segment mapping matrix, laying the data foundation for subsequent pruning strategies. However, when processing samples containing extremely long read segments, the static threshold may cause the sliding window pruning to be over-triggered, resulting in indiscriminate truncation of functional segments, coverage gaps during the alignment stage, and consequently, missing variant fields, affecting the accuracy of population diversity analysis. Therefore, step S3 needs to extract a sensitive region index set at risk of over-pruning based on the mapping matrix by analyzing the triggering frequency, segment triggering dispersion, and coverage gap balance of the sliding window pruning, in order to support adaptive adjustment of the pruning strategy.

[0073] Step S3 analyzes the triggering characteristics of sliding window pruning in different read length intervals based on the high and low quality read segment mapping matrix generated in step S2, calculates the segment triggering dispersion and coverage gap balance, constructs a pruning risk mapping model, and finally extracts the sensitive area index set with over-truncation risk.

[0074] Step S3 includes the following:

[0075] S3.1, Statistical sliding window trigger frequency sequence:

[0076] Within each read length segment, a sliding window pruning process is simulated, recording the number of times the pruning operation is triggered and the total number of operations. The trigger frequency is defined as the ratio of the number of pruning triggers to the total number of operations, reflecting the frequency of the pruning operation within the corresponding read length segment. After statistics are completed, a sequence containing the trigger frequency for each read length segment is generated.

[0077] The statistical sliding window trigger frequency sequence aims to quantify the execution of pruning operations across different read length segments, facilitating the analysis of the distribution characteristics of pruning operations. By simulating the pruning process and recording triggering events, the frequency differences of pruning operations can be accurately captured, thereby identifying read length segments with abnormal frequencies and improving the targeting and accuracy of the analysis.

[0078] S3.2, Calculate the segment trigger dispersion:

[0079] A change-point potential energy analysis is performed on the trigger frequency sequence to identify change points and divide the sequence into multiple sub-intervals. Within each sub-interval, the difference between the maximum and minimum trigger frequency values ​​is calculated to obtain the frequency fluctuation amplitude. Then, the frequency fluctuation amplitude of each sub-interval is multiplied by the sub-interval length, and the results for all sub-intervals are summed to obtain the segment trigger dispersion. The segment trigger dispersion measures the degree of fluctuation of the trigger frequency between different read length segments, reflecting the non-uniformity of clipped triggering.

[0080] Calculating the segment trigger dispersion is used to assess the fluctuation of the trigger frequency and identify read length segments with uneven clipping operations. By employing variable-point potential energy analysis combined with weighted fluctuation amplitude calculation, local variation characteristics of the trigger frequency sequence can be captured, facilitating the identification of regions with abnormal clipping triggers, thereby improving the accuracy and reliability of the analysis.

[0081] S3.3, Calculate the coverage gap balance:

[0082] The coverage surface of each reading segment is extracted from the comparison curves before and after trimming, and the peak height and valley depth in the coverage surface are analyzed. The coverage gap balance is defined as 1 minus the ratio of the difference between the peak height and valley depth to the total height of the coverage surface (the sum of the peak height and valley depth), which is used to measure the gap symmetry of the coverage surface and reflect the impact of the trimming operation on the coverage integrity.

[0083] Calculating the coverage gap evenness aims to assess the impact of trimming operations on the integrity of the coverage surface and identify areas of uneven coverage caused by trimming. By analyzing the difference between peak height and valley depth, the degree of change in the coverage surface can be quantified, facilitating the assessment of the rationality of the trimming operation and thus improving the accuracy and consistency of data processing. The processing focuses on the extraction and difference analysis of surface features to ensure that the evenness truly reflects the trimming effect and supports subsequent risk assessment.

[0084] S3.4, Construct a trimming risk mapping model:

[0085] A two-dimensional coordinate system is constructed using the segment trigger dispersion and coverage gap balance as orthogonal axes, marking the coordinate points of each read length segment. A two-dimensional Laplace transform is applied to the coordinate point set to generate a potential energy distribution surface. On the potential energy distribution surface, the minimum closed loop is extracted along the zero potential line, and the area enclosed by the minimum closed loop is calculated, defined as the truncation sensitivity coefficient. The magnitude of the truncation sensitivity coefficient reflects the level of over-truncation risk; a larger area indicates a higher risk.

[0086] The purpose of constructing a pruning risk mapping model is to integrate segment trigger dispersion and coverage gap balance into a single risk indicator to comprehensively assess the risk of pruning operations. By using two-dimensional Laplace transform and closed-loop area calculation, the influence of dispersion and balance can be combined to quantify the degree of risk and facilitate the identification of high-risk read length segments.

[0087] S3.5, Extract the sensitive region index set:

[0088] For each read length segment, if its truncation sensitivity coefficient is higher than a set threshold, it is included in the sensitive region index set. The sensitive region index set contains read length segments that are at risk of excessive truncation, and is used to adjust subsequent pruning strategies.

[0089] Extracting the sensitive region index set aims to identify read length segments at risk of excessive truncation, providing a basis for optimizing pruning strategies. Setting thresholds to filter high-risk segments ensures the relevance and effectiveness of the analysis, reducing the negative impact of pruning operations on the data.

[0090] Step S3 achieves accurate identification of over-truncation risk read length segments by statistically analyzing the sliding window trigger frequency sequence, calculating the segment trigger dispersion and coverage gap balance, constructing a pruning risk mapping model, and extracting a sensitive region index set. The entire process addresses the characteristics of ultra-long read fragment samples, resolving the issues of coverage gaps and missing variant fields caused by static thresholds, and providing technical support for the scientific adjustment of pruning strategies.

[0091] Step S3 analyzes the trigger frequency, segment trigger dispersion, and coverage gap balance of the sliding window pruning based on the mapping matrix, extracting a set of sensitive regions at risk of over-truncation. These steps lay the foundation for identifying high-risk segments. However, when processing samples containing excessively long read segments, static thresholds may lead to indiscriminate truncation of functional segments, affecting the integrity of alignment coverage and variant detection. Therefore, step S4 requires adaptive pruning for the sensitive region index set through backtracking and adjusting the pruning strategy to reduce over-truncation and maintain data integrity.

[0092] Step S4 involves backtracking the sensitive region index set to the original read length segment position, generating a pruning and release adjustment table, injecting it into the pruning task configuration, and replacing the static threshold.

[0093] Step S4 includes the following:

[0094] S4.1, backtrack the sensitive region index set to the original read length segment position:

[0095] The sensitive region index set contains several sensitive read segments, mapping each segment to its corresponding location in the original sequencing file. Each sensitive read segment is located by its start and end positions in the original sequencing file. Simultaneously, by combining quantile tags, the sequence beginning and end features of each sensitive read segment are associated to generate a sensitive region location table. The sensitive region location table records the start and end positions and quantile tags of each sensitive read segment.

[0096] The purpose of tracing the sensitive region index back to the original read segment location is to provide precise location data for subsequent pruning adjustments, ensuring that pruning operations can be accurately applied to high-risk regions. By establishing the association between sensitive read segments and specific locations in the original sequencing file, the accuracy and consistency of data processing can be improved, avoiding erroneous operations on non-target segments and thus ensuring the reliability of the processing.

[0097] S4.2, Generate the trimming and slow-release adjustment table:

[0098] For each sensitive read length segment, dynamic pruning parameters are calculated to reduce the risk of over-truncation. Using the quality stability score and the truncation sensitivity coefficient, the sum of the quality stability score and the truncation sensitivity coefficient plus 1 yields a release factor used to adjust the pruning strictness. The quality stability score reflects the quality consistency of the read length segment; a higher value indicates greater quality stability. The truncation sensitivity coefficient represents the degree of risk of over-truncation; a higher value indicates greater risk. By using the quality stability score as the numerator and the truncation sensitivity coefficient plus 1 as the denominator, the release factor balances the relationship between quality stability and truncation risk; a larger value indicates a greater need to relax the pruning strictness to reduce over-truncation of high-quality segments.

[0099] Next, for each sensitive read segment, its release factor and original static threshold are obtained. The original static threshold is multiplied by the release factor plus 1 to obtain the dynamically adjusted threshold, and the result is recorded in the trimmed release adjustment table. The trimmed release adjustment table includes the identifier of the sensitive read segment, its start position, its end position, and its corresponding dynamic threshold.

[0100] The purpose of generating a pruning buffer adjustment table is to provide personalized pruning parameters for sensitive read length segments, thereby reducing data loss caused by uniform truncation. By introducing buffer factors and dynamic thresholds, the pruning strategy can be adjusted according to the quality characteristics and risk level of the segment, enhancing the flexibility and adaptability of processing and improving the integrity of data retention.

[0101] S4.3, Inject trimming task configuration:

[0102] The dynamic thresholds in the trimming buffer adjustment table are integrated with the original static thresholds to generate a new trimming task configuration. In this new configuration, dynamic thresholds are applied to sensitive read length segments, while the original static thresholds are applied to non-sensitive read length segments. Subsequently, the trimming task configuration file is updated to ensure that the trimming tool can identify and apply the appropriate thresholds based on the start and end positions of the segments during execution.

[0103] The purpose of injecting pruning task configuration is to enable adaptive adjustment of pruning strategies, allowing pruning operations to be differentiated based on the risk level of segments. By integrating dynamic thresholds and original static thresholds, the targeting and effectiveness of pruning operations can be enhanced, reducing unnecessary data truncation caused by uniform thresholds.

[0104] S4.4, Replace the static threshold:

[0105] Inject the updated trimming task configuration into the trimming tool's configuration file to ensure that the trimming tool uses dynamic thresholds for sensitive read length segments and the original static thresholds for non-sensitive read length segments during runtime. After updating the trimming task configuration, the trimming tool can then perform processing according to the adaptive trimming strategy.

[0106] The purpose of replacing the static threshold is to enable the pruning tool to adaptively perform pruning operations based on the risk level of the segment. By injecting the updated configuration into the pruning tool, the flexibility and accuracy of processing can be improved, the impact of over-truncation on sensitive read-length segments can be reduced, and thus the stability of data processing can be enhanced.

[0107] Step S4 achieves adaptive adjustment of the pruning strategy by sequentially executing four sub-steps: backtracking the sensitive region index set to the original read length segment position, generating the pruning slow-release adjustment table, injecting the pruning task configuration, and replacing the static threshold. For samples containing ultra-long read fragments, this process effectively reduces indiscriminate truncation of functional segments, improves the accuracy of pruning operations and the integrity of data preservation, and provides technical support for the stability and accuracy of whole-genome genetic diversity cloud analysis.

[0108] Steps S1 to S4, described above, addressed the issue of indiscriminate truncation of functional regions in ultra-long read samples due to static thresholds by collecting raw sequencing files, dividing read length segments, constructing a quality mapping matrix, identifying sensitive region index sets, and generating a pruning and release adjustment table. These steps also updated the pruning task configuration. However, simply completing the pruning adjustment is insufficient to guarantee data integrity and reliability. Further verification of the coverage stability and variant field integrity of the pruned sequences in sensitive regions is necessary to avoid coverage gaps during alignment and field omissions in variant detection, ensuring the accuracy of subsequent population diversity analysis. Therefore, step S5 focuses on performing the adjusted pruning and, through alignment and verification processes, quantitatively evaluating the data quality of sensitive regions to refine the adaptive pruning strategy.

[0109] Step S5, based on the pruning task configuration output in step S4, performs the adjusted pruning and compares the output sequence, extracts the coverage curve and mutation distribution, verifies the coverage stability of sensitive segments and the integrity of the mutation field, and finally writes the verification results into the parameter cache table.

[0110] Step S5 includes the following:

[0111] S5.1, Perform the adjusted cropping:

[0112] The Trimmomatic tool was used to trim the original sequencing files according to the trimming task configuration. The trimming task configuration included adjustment parameters for dynamic and static thresholds. The dynamic threshold was adaptively adjusted based on the characteristics of sensitive read regions, while the static threshold served as a fixed standard for non-sensitive read regions. After trimming, the trimmed sequencing files were generated.

[0113] S5.2, Compare the cropped sequence:

[0114] The alignment tool BWA was used to align the trimmed sequencing file with the reference genome, generating an alignment file. This file records the alignment position and coverage of the trimmed sequence on the reference genome, fully reflecting the matching relationship between the sequence and the reference genome.

[0115] S5.3, Extract the coverage curve:

[0116] The coverage depth sequence of sensitive segments is extracted from the alignment file, and the coverage curve fluctuation coefficient of each sensitive segment is calculated. The calculation method of the coverage curve fluctuation coefficient is as follows: first, the maximum and minimum values ​​of the coverage depth sequence are determined, the difference between the two is calculated, and then the difference is divided by the mean of the coverage depth sequence to obtain the coverage curve fluctuation coefficient, which is used to quantify the degree of fluctuation of the coverage curve.

[0117] The goal of extracting coverage curves is to assess the coverage distribution of sensitive regions on the reference genome after pruning and to identify areas of uneven coverage. By calculating the coverage curve fluctuation coefficient, the stability of coverage can be quantitatively described, facilitating the determination of whether the pruning operation resulted in coverage gaps, thereby ensuring the reliability of data analysis.

[0118] S5.4, verify coverage stability:

[0119] The coverage curve fluctuation coefficient of each sensitive section is compared with a preset stability threshold to determine coverage stability. If the coverage curve fluctuation coefficient is less than the preset stability threshold, the coverage stability of the sensitive section is determined to meet the requirements; if the coverage curve fluctuation coefficient is greater than or equal to the preset stability threshold, the coverage stability is determined to be insufficient.

[0120] The purpose of verifying coverage stability is to ensure that the coverage curve of sensitive segments after pruning does not exhibit abnormal fluctuations, thus avoiding the impact of uneven coverage on the accuracy of subsequent analyses. By comparing with a preset stability threshold, unstable coverage segments can be quickly identified, thereby improving the reliability of data quality.

[0121] S5.5, Extract the variation distribution:

[0122] The file was compared using the mutation detection tool GATK to perform mutation detection, and mutation detection results were generated. For sensitive sections, the completeness of the mutation fields in the mutation detection results was checked, specifically including whether there were any omissions in the mutation type, location, and related parameters, and the check results were recorded.

[0123] The purpose of extracting the variant distribution is to evaluate the quality of variant detection in sensitive regions after pruning and to ensure the integrity of variant information. By checking the integrity of the variant field, we can identify variant information loss caused by the pruning operation, thereby improving the accuracy and reliability of variant detection.

[0124] S5.6, Verify the integrity of the variant fields:

[0125] Based on the integrity check results of the variant fields, determine whether the variant fields of the sensitive section are complete. If there are no missing variant fields in the variant detection results, the variant fields are considered complete; if there are missing variant types, locations, or related parameters, the variant fields are considered incomplete.

[0126] The goal of verifying the integrity of the variation field is to ensure that variation information in sensitive regions is not lost after pruning, thus avoiding the impact of missing fields on the accuracy of population genetic diversity analysis. Directly checking the integrity of the variation field can improve the accuracy and consistency of data processing.

[0127] S5.7, Write to the parameter cache table:

[0128] The coverage stability and variation integrity status of each sensitive segment are recorded in the parameter cache table. Based on the verification results, if the coverage curve fluctuation coefficient of a sensitive segment is less than the preset stability threshold and the variation field is complete, it is marked as "verification passed"; if the coverage curve fluctuation coefficient is greater than or equal to the preset stability threshold, or the variation field is missing, it is marked as "adjustment required".

[0129] Step S5, by performing adjusted pruning, sequence alignment, extraction of coverage curves and variant distributions, and verification of coverage stability and variant field integrity in sensitive regions, achieves a comprehensive assessment of the quality of the pruned data. For the characteristics of ultra-long read samples, this step ensures the reliability of the pruning results in sequence alignment and variant detection. Finally, the verification results are recorded through a parameter cache table, providing technical support for the stability and accuracy of whole-genome genetic diversity cloud analysis.

[0130] The above formulas are all dimensionless calculations. The formulas are derived from software simulations based on a large amount of collected data to obtain the most recent real-world results. The preset parameters in the formulas are set by those skilled in the art according to the actual situation.

[0131] It should be noted that the system of the present invention can be deployed on the device itself to realize embedded applications, or it can run on a PC or other terminal with a user interface, thereby meeting a variety of hardware environments and usage requirements.

[0132] The foregoing has only described certain exemplary embodiments of the present invention by way of illustration. Undoubtedly, those skilled in the art can modify the described embodiments in various ways without departing from the spirit and scope of the present invention. Therefore, the foregoing drawings and descriptions are illustrative in nature and should not be construed as limiting the scope of protection of the claims of the present invention.

[0133] It should be noted that, in this document, the use of relational terms such as "first" and "second" is merely to distinguish one entity or operation from another, and does not necessarily require or imply any such actual relationship or order between these entities or operations. Furthermore, the terms "comprising," "including," or any other variations thereof are intended to cover non-exclusive inclusion, such that a process, method, article, or apparatus that comprises a list of elements includes not only those elements but also other elements not expressly listed, or elements inherent to such a process, method, article, or apparatus. Without further limitations, an element defined by the phrase "comprising one..." does not exclude the presence of other identical elements in the process, method, article, or apparatus that includes the element.

[0134] The above description is merely a specific embodiment of this application, but the scope of protection of this application is not limited thereto. Any variations or substitutions that can be easily conceived by those skilled in the art within the scope of the technology disclosed in this application should be included within the scope of protection of this application. Therefore, the scope of protection of this application should be determined by the scope of the claims.

Claims

1. A cloud-based method for whole-genome genetic diversity analysis based on workflow automation, characterized in that: Including the following steps: S1: By analyzing the read length distribution and the first and last base quality characteristics of the original sequencing files, the read lengths are divided into different segments and corresponding quantile tags are generated; S2: By calling the quality histogram distribution output by FastQC and combining it with quantile labels, a high- and low-quality read segment mapping matrix is ​​constructed to mark the quality stability of each read length segment; The raw sequencing files were analyzed by FastQC to generate a quality histogram distribution report, quality characteristic indicators were extracted, a high- and low-quality read mapping matrix was constructed using quantile tags, and stability scores were calculated by setting thresholds and nonlinear fusion to mark the quality stability of each read segment. Finally, the mapping matrix and stability tags were stored as structured data. S3: Determine the trigger frequency of sliding window pruning in different read length intervals based on the mapping matrix, calculate the segment trigger dispersion and coverage gap balance, synthesize the truncation sensitivity coefficient through the preset pruning risk mapping model, and extract the sensitive area index set; By analyzing the mapping matrix, the trigger frequency of sliding window pruning in different read length intervals is determined, the segment trigger dispersion and coverage gap balance are calculated, and the truncation sensitivity coefficient is synthesized using the preset pruning risk mapping model. Read length segments with truncation sensitivity coefficients higher than the set threshold are extracted to form a sensitive area index set. S4: Backtrack the sensitive region index set to the original read length segment position, generate a pruning and release adjustment table and inject the pruning task configuration to replace the static threshold; A sensitive region location table is generated by tracing back the sensitive region index set to the start and end positions of sensitive reads in the original sequencing file; a sustained release factor is calculated using the quality stability score and truncation sensitivity coefficient, and a dynamic threshold is generated by combining the original static threshold to construct a pruning and sustained release adjustment table; Integrate the dynamic thresholds and original static thresholds in the clipping release adjustment table into the clipping task configuration, update the clipping tool configuration file to apply dynamic thresholds to sensitive read length segments and apply original static thresholds to non-sensitive read length segments; S5: Perform the adjusted pruning and compare the output sequence, extract the coverage curve and variation distribution, verify the coverage stability of sensitive segments and the integrity of the variation field, and then write the results to the parameter cache table.

2. The cloud-based method for whole-genome genetic diversity analysis based on workflow automation as described in claim 1, characterized in that, Step S1 includes the following: Raw sequencing files in FASTQ format were collected, and read length distribution was analyzed to determine minimum length, maximum length, and distribution characteristics. Reads were segmented into intervals based on the length and quality characteristics of the start and end bases, and quantile labels were assigned to each interval according to length quantiles and quality categories.

3. The cloud-based method for whole-genome genetic diversity analysis based on workflow automation according to claim 2, characterized in that, Step S2 includes the following: Quality characteristics include average quality score, quality score fluctuation range, GC content deviation, and repetitive sequence ratio; the average quality score is obtained by calculating the average of the quality scores of all bases within the read length segment; The fluctuation range of the quality score is obtained by calculating the difference between the maximum and minimum quality scores within the read length segment; the deviation of GC content is determined by calculating the absolute percentage deviation between the GC content of the read length segment and the global average GC content; the proportion of repeated sequences is obtained by calculating the proportion of repeated sequences within the read length segment to the total sequence.

4. The cloud-based method for whole-genome genetic diversity analysis based on workflow automation according to claim 3, characterized in that, Step S3 includes the following: The processing logic for segment trigger dispersion is as follows: the distribution of trigger events is obtained by using the sliding window trigger frequency sequence of each read length segment in the statistical mapping matrix. Then, the frequency sequence is subjected to variable point potential energy analysis to identify key points of frequency change. The dispersion index is obtained by calculating the fluctuation variance, thereby quantitatively describing the uneven distribution of clipping triggers.

5. The cloud-based method for whole-genome genetic diversity analysis based on workflow automation according to claim 4, characterized in that, Step S3 also Includes the following: The processing logic for coverage gap balance is as follows: extract the segment coverage surface from the comparison curves before and after trimming, analyze the change characteristics of the coverage range, and then perform a peak-valley symmetry measurement on the segment coverage surface. By calculating the degree of symmetry between the peaks and valleys, the gap balance index is obtained, thereby assessing the impact of the trimming operation on the coverage integrity.

6. The cloud-based method for whole-genome genetic diversity analysis based on workflow automation according to claim 5, characterized in that, Step S5 Includes the following: The original sequencing file is trimmed using the dynamic and static threshold adjustment parameters configured in the trimming task, generating a trimmed sequencing file. The trimmed sequencing file is then aligned with the reference genome using the alignment tool BWA, generating an alignment file. Coverage depth sequences of sensitive regions are extracted from the alignment file, and the coverage curve fluctuation coefficient is calculated to quantify the degree of fluctuation of the coverage curve. The coverage curve fluctuation coefficient is compared with a preset stability threshold to determine the coverage stability of sensitive regions.

7. The cloud-based method for whole-genome genetic diversity analysis based on workflow automation according to claim 6, characterized in that, Step S5 also Includes the following: The GATK mutation detection tool is used to compare and detect mutations in the file, generate mutation detection results, and check the integrity of the mutation fields in sensitive sections. Based on the integrity check results of the mutation fields, the integrity of the mutation fields in sensitive sections is determined. The coverage stability and mutation integrity status of sensitive sections are recorded in the parameter cache table, and the status is marked according to the verification results.

Images