A joint alignment and quantification method for large-scale mass spectrometry cohort
By employing a joint alignment and quantification method using end-to-end deep learning networks, the problems of error accumulation and lack of coordination in large-scale mass spectrometry queues are solved, enabling efficient and accurate data processing and analysis.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- SHANGHAI DEV CENT OF COMP SOFTWARE TECH
- Filing Date
- 2026-04-24
- Publication Date
- 2026-07-03
AI Technical Summary
Existing high-throughput mass spectrometry data processing workflows suffer from error accumulation and lack of coordination when handling large-scale queues, including defects such as error propagation, lag in intensity drift correction, and inconsistency in objective functions.
We employ a unified backpropagation process using an end-to-end deep learning network. By simultaneously optimizing retention time alignment, cross-batch intensity drift correction, and quantitative regression through a joint loss function, we achieve joint alignment and quantification of data using a one-dimensional monotonic radial basis function space transformation network, an explicit intensity drift correction layer, and a two-stream attention mechanism.
It achieves efficient and accurate alignment and quantification of large-scale mass spectrometry cohort data, overcomes the shortcomings of traditional step-by-step strategies, reduces data processing time and storage space consumption, and improves the accuracy and consistency of analysis.
Smart Images

Figure CN122112596B_ABST
Abstract
Description
Technical Field
[0001] This application relates to the interdisciplinary field of bioinformatics, artificial intelligence and high-performance computing, and in particular to a joint alignment and quantification method for large-scale mass spectrometry cohorts. Background Technology
[0002] Chromatography-mass spectrometry (GC-MS) techniques, such as LC-MS and GC-MS, have become core tools for proteomics, metabolomics, lipidomics, and clinical drug analysis due to their high throughput, high sensitivity, and wide dynamic range. With the rise of precision medicine and large-scale cohort studies (such as serum / urine sample cohorts of thousands of cases), mass spectrometry data analysis is facing a paradigm shift from "small-sample fine-grained analysis" to "large-scale cohort industrial analysis." However, existing high-throughput mass spectrometry data processing workflows (such as XCMS, MaxQuant, and OpenMS) suffer from at least error accumulation and lack of coordination when processing large-scale cohorts due to their "multi-step serial" nature. Specifically, existing processing workflows are typically designed as disjointed linear steps: first, retention time (RT) alignment; second, peak detection and matching; and finally, intensity normalization and quantification. This step-by-step strategy has significant drawbacks, including error propagation, lag in intensity drift correction, and inconsistent objective functions. Summary of the Invention
[0003] The purpose of this application is to provide a joint alignment and quantification method for large-scale mass spectrometry cohorts, which can overcome the problems of error accumulation and lack of coordination in existing technologies.
[0004] To achieve the above objectives, this application provides the following solution:
[0005] Firstly, this application provides a joint alignment and quantification method for large-scale mass spectrometry queues, which outputs a quantitative feature matrix through an end-to-end deep learning network; wherein, during the unified backpropagation process of the end-to-end deep learning network, retention time alignment, cross-batch intensity drift correction, and quantitative regression are simultaneously optimized based on a joint loss function, and the training steps of the end-to-end deep learning network include:
[0006] The original three-dimensional feature tensors of all batches are obtained. The original three-dimensional feature tensors of one batch are used as standard samples, and the original three-dimensional feature tensors of the remaining batches are aligned with the standard samples. The original three-dimensional feature tensors are extracted from the chromatography-mass spectrometry queue.
[0007] A one-dimensional monotonic radial basis function space transformation network is used to learn a time transformation function from aligning the time grid to preserving the time axis for each sample to achieve time alignment of the original three-dimensional feature tensors. Intensity compensation is performed on the original three-dimensional feature tensors based on the time transformation function and its first derivative. The intensity compensation is used to maintain the conservation of the chromatographic peak integral area when the time axis is nonlinearly stretched. The original three-dimensional feature tensors after time alignment and intensity compensation are intermediate three-dimensional feature tensors.
[0008] Based on the sample batch identifier, injection order, and acquisition time of the intermediate three-dimensional feature tensor, a scaling factor and a bias factor are generated through an explicit intensity drift correction layer. An element-wise affine transformation is then performed on the intermediate three-dimensional feature tensor based on the scaling factor and the bias factor to achieve intensity drift correction, thereby obtaining the target three-dimensional feature tensor.
[0009] The target 3D feature tensor is enhanced by a dual-stream attention mechanism based on global flow features and differential flow features. The enhanced target 3D feature tensor is then input into a quantitative regression module based on sequence modeling, which outputs the quantitative feature matrix.
[0010] Secondly, this application provides a computer device, including: a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the aforementioned joint alignment and quantification method for large-scale mass spectrometry queues.
[0011] Thirdly, this application provides a computer program product, including a computer program that, when executed by a processor, implements the aforementioned joint alignment and quantification method for large-scale mass spectrometry queues.
[0012] According to the specific embodiments provided in this application, the following technical effects are disclosed:
[0013] This application provides a joint alignment and quantification method for large-scale mass spectrometry cohorts. It employs an end-to-end deep learning network to analyze data from large-scale chromatography-mass spectrometry cohorts. During the unified backpropagation process of the end-to-end deep learning network, retention time alignment, cross-batch intensity drift correction, and quantification regression are performed simultaneously. This means that retention time alignment, cross-batch intensity drift correction, and quantification regression are treated as a unified joint operation, rather than sequential operations, overcoming the shortcomings of step-by-step strategies in traditional techniques. Specifically, retention time alignment, cross-batch intensity drift correction, and quantification regression are implemented through different network modules within the end-to-end deep learning network. Retention time alignment is achieved through a one-dimensional monotonic radial basis function space transformation network; cross-batch intensity drift correction is achieved through an explicit differentiable intensity drift correction layer; and quantification regression is achieved through a two-stream attention mechanism and a quantification regression module. The network parameters of these modules are all synchronously optimized based on a joint loss function. Attached Figure Description
[0014] To more clearly illustrate the technical solutions in the embodiments of this application or the prior art, the drawings used in the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of this application. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.
[0015] Figure 1 This is a schematic diagram of the overall process of the joint alignment and quantification method for large-scale mass spectrometry queues provided in the embodiments of the present invention.
[0016] Figure 2 This is a detailed architectural block diagram of the end-to-end deep learning network in an embodiment of the present invention.
[0017] Figure 3 This is a schematic diagram illustrating the data mapping and end-to-end sparse computing principle of the virtualized columnar interface in this embodiment of the invention.
[0018] Figure 4 This is a schematic diagram of the differentiable time transformation and peak area conservation compensation principle based on one-dimensional RBF in an embodiment of the present invention.
[0019] Figure 5 This is a schematic diagram illustrating the explicit intensity drift correction and QC manifold alignment principle in an embodiment of the present invention.
[0020] Figure 6 This is a schematic diagram of the boundary detection and integral quantification principle based on hybrid supervision in an embodiment of the present invention. Detailed Implementation
[0021] The technical solutions of the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of this application, and not all embodiments. Based on the embodiments of this application, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of this application.
[0022] To make the above-mentioned objectives, features and advantages of this application more apparent and understandable, the application will be further described in detail below with reference to the accompanying drawings and specific embodiments.
[0023] The joint alignment and quantification method for large-scale mass spectrometry cohorts described in this application can be applied to terminal device servers. The terminal device can be a host computer, laptop computer, or other computer equipment, and the server can be a single server or a server cluster. The following describes a specific embodiment of the joint alignment and quantification method for large-scale mass spectrometry cohorts, using an application to a terminal device as an example.
[0024] One embodiment of this application provides an adaptive alignment and joint quantification method for large-scale chromatography-mass spectrometry cohorts (including but not limited to liquid chromatography-mass spectrometry (LC-MS) and gas chromatography-mass spectrometry (GC-MS)). The overall process is described in [link to documentation]. Figure 1 The heterogeneous raw data pool is processed, and a quantitative feature matrix (i.e., a standardized quantitative matrix) is output through an end-to-end deep learning network; where the end-to-end deep learning network (corresponding to...) Figure 1 The deep learning joint optimization network (in the context of deep learning) simultaneously optimizes retention time alignment, cross-batch intensity drift correction, and quantitative regression based on a joint loss function during the unified backpropagation process. The training steps of the end-to-end deep learning network include steps S110, S120, S130, and S140. The detailed architecture diagram of the end-to-end deep learning network is shown below. Figure 2 .
[0025] Step S110: Obtain the original three-dimensional feature tensors of all batches, use the original three-dimensional feature tensors of one batch as the standard sample, and align the remaining batches with the sample of that batch as the standard (i.e., align the original three-dimensional feature tensors of the remaining batches with the standard sample as the standard). The original three-dimensional feature tensors are extracted from the chromatography-mass spectrometry queue.
[0026] In this preferred step, obtaining the original three-dimensional feature tensors of all batches specifically includes: parsing the metadata index of the original mass spectrometry data file, and establishing an index mapping relationship from physical storage address to logical feature space based on the metadata index (corresponding to...). Figure 1In step S1 of the virtualization index construction, without decoding all the raw mass spectrometry data files, data blocks are read on demand from heterogeneous row-oriented or column-oriented storage files on the disk through index mapping relationships, and multiple data blocks are assembled into a batch of raw three-dimensional feature tensors in a streaming manner (corresponding to...). Figure 1 (Streaming sparse micro-batch loading in step S2).
[0027] For example, the raw mass spectrometry data file includes: (a) a raw mass spectrometry data file stored in row format and / or (b) a data file in column format or block format obtained by preprocessing the raw mass spectrometry data file.
[0028] The virtualized columnar data access interface is configured as follows: when the input is a raw mass spectrometry data file stored in row format, an external index is generated based on the scan / spectrum-level index structure, and data blocks are located and read on demand (row-based on-demand loading mode); when the input is a data file in columnar or block storage format, the corresponding columnar data access interface is called and its metadata index is used to locate the target data block (column-based block reading mode); without decoding all spectrum data, data blocks in heterogeneous row-based storage files on disk are read on demand through index mapping, and multi-sample data are assembled in a streaming manner into a micro-batch three-dimensional feature tensor or its sparse equivalent representation as a batch of raw three-dimensional feature tensors.
[0029] The virtualized columnar data access interface selects or switches between row-based on-demand loading mode and columnar block reading mode based on at least one of the following: the storage format identifier of the input file, the training round number threshold, the storage medium bandwidth, or the cache hit rate.
[0030] The original three-dimensional feature tensor is denoted as ,in, For sample batch dimension, To retain the time dimension for discrete elements, The feature dimension refers to the dimensions of the original three-dimensional feature tensor, which include the sample batch dimension, the discrete retention time dimension, and the feature dimension. The feature dimension is defined as follows: dividing the preset mass-to-charge ratio range into multiple non-overlapping or partially overlapping sub-intervals, with each sub-interval corresponding to one or more indices in the feature dimension; or mapping the pre-extracted spectral features to a set of channels on the feature dimension. The pre-extracted spectral features include pre-identified regions of interest, specific isotope distribution clusters, and fragment ion groups.
[0031] Elements in constructing the original three-dimensional feature tensor At that time, the data block is placed in the time window Inside and fall into the feature unit fThe polymerization intensity value is obtained by polymerizing one or more ion signal intensities within the corresponding mass-to-charge ratio range, and the element is determined based on the polymerization intensity value. .
[0032] During the aggregation process, the aggregation operators used include at least one of summation, maximum value taking, weighted average, and peak area integration.
[0033] Determining elements based on polymerization strength value Specifically, this includes: performing logarithmic transformation, range normalization, or standard deviation standardization on the aggregate intensity values to improve the comparability of different original 3D feature tensors in terms of feature dimensions and suppress the impact of extreme intensity values on training; and using the processed aggregate intensity values as elements. .
[0034] For example, without decoding the entire raw mass spectrometry data file, data blocks are read on demand from a heterogeneous row-oriented storage file on disk using an index mapping relationship, and multiple data blocks are assembled into a batch of raw three-dimensional feature tensors in a streaming manner, specifically including:
[0035] By calling the vendor's raw data access library or a general parsing adapter, the index list, offset table, scan directory tree, or equivalent structure in the raw mass spectrometry data file is parsed. The retention time, scan number, total ion current summary, and byte start offset and byte length in the physical file corresponding to each scan are extracted from the structure without decoding the complete m / z-intensity array. Based on the extracted metadata, a memory-resident lookup table and / or a persistent index file are constructed, mapping sample identifiers and retention time window identifiers to physical file paths and data block locations. During the training or inference phase, the file location offset of the required data blocks is calculated according to the index table. A lazy loading method is adopted, reading only the part of the data required for the current micro-batch training or inference from the disk, thereby avoiding loading the entire file into memory.
[0036] Optionally, data blocks can be sorted by file offset and adjacent or overlapping data blocks can be merged into consecutive read requests to reduce random I / O.
[0037] Optionally, an asynchronous prefetching and double buffering mechanism can be used to preload the data blocks required for the next microbatch while processing the current microbatch.
[0038] When the raw mass spectrometry data file is in a pre-converted columnar or block storage format such as HDF5, Zarr, Lance, or Parquet, the virtualized columnar data access interface is configured to: directly utilize the built-in columnar access and block compression mechanisms of the aforementioned formats for data reading; the chunk / block layout of the columnar or block storage format data file is organized according to sample identifier, retention time window, and / or mass-to-charge ratio bin (or ROI feature channel), ensuring that a single read covers the minimum set of contiguous data blocks required for the current microbatch. The target sample, retention time period, and mass-to-charge ratio bin are determined using metadata or a built-in index; multiple data blocks are loaded and decompressed in parallel using a multi-threaded prefetch queue; and the prefetch queue length and buffer size are dynamically adjusted based on hardware resource status to maximize I / O throughput and reduce access latency.
[0039] Furthermore, step S110 also includes: converting the raw mass spectrometry data file stored in row format into a columnar or block-based storage format data file, and generating metadata and index information consistent with the logical feature space during the conversion process, which are used for subsequent location and reading of data blocks according to retention time windows and / or mass-to-charge ratio bins.
[0040] The joint alignment and quantification method can employ a memory-GPU memory hierarchical management and end-to-end sparse computation strategy. This involves reading data blocks from disk on demand through index mapping relationships, assembling multi-sample data into micro-batch 3D feature tensors in a streaming manner, and illustrating the data mapping and end-to-end sparse computation principle using a virtualized columnar interface as shown below. Figure 3 As shown. This strategy specifically includes:
[0041] (1) Streaming micro-batch loading (I / O level optimization): The retention time sliding window mechanism and sample micro-batch strategy are adopted. Only a portion of the sample data that is within the retention time sliding window required for the current calculation is loaded into memory, instead of loading the complete sample data all at once. The processing of the complete sample data is completed through loop iteration.
[0042] (2) Sparse compression and dynamic mapping (storage-level optimization): Micro-batch tensors are stored in coordinate list (COO) or compressed sparse row (CSR) format. Through dynamic ROI mapping technology, hash mapping or index lookup table is used to map high-dimensional sparse physical m / z coordinates into compact logical feature channel indexes, thereby filtering out more than 99% of invalid background zero values before entering the computation graph.
[0043] (3) Sparsity-Aware Computation (Operator-Level Optimization): Sparsity-Aware Computation is performed in the network module to construct a binary mask corresponding to the non-zero features. The mask is used to shield invalid regions in the Softmax step of the Transformer self-attention mechanism, and the PaddingMask mechanism is used to block gradient updates of zero-value regions during backpropagation, so as to focus on the learning of chemical features of effective ions.
[0044] (4) Adaptive resource scheduling (system-level optimization): Real-time monitoring of GPU memory usage. When memory usage exceeds a preset threshold, a dynamic adjustment mechanism is automatically triggered to perform at least one of the following operations: reducing the micro-batch size, shortening the time window length, or dynamically increasing the sparsity noise threshold, in order to prevent memory overflow (OOM).
[0045] Step S120: A time transformation function from aligning the time grid to preserving the time axis is learned for each original three-dimensional feature tensor through a one-dimensional monotonic radial basis function space transformation network to achieve time alignment of each original three-dimensional feature tensor. Intensity compensation is performed on the original three-dimensional feature tensor based on the time transformation function and its first derivative. The intensity compensation is used to maintain the conservation of the integral area of the chromatographic peak when the time axis is nonlinearly stretched. The original three-dimensional feature tensor after time alignment and intensity compensation is the intermediate three-dimensional feature tensor.
[0046] Among them, the one-dimensional monotonic radial basis function space transformation network includes a localization subnetwork (i.e. Figure 2 The system comprises a localization subnetwork, a mesh generator, and a differentiable sampler. The localization subnetwork, based on the input original 3D feature tensor, outputs transformation parameters characterizing the time transformation function. These transformation parameters include a set of time control points defined in the normalized time domain [0,1]. The positions of these time control points are adaptively determined based on the peak density and / or gradient program variation nodes and / or prior anchor points of the total ion current curve. The time control points are either fixed parameters or finely adjusted within a preset range using learnable offsets while maintaining monotonicity. The mesh generator generates mapped coordinates between the aligned time mesh and the original retained time axis based on the transformation parameters. The differentiable sampler performs differentiable interpolation sampling on the original 3D feature tensor in the retained time dimension based on the mapped coordinates to obtain geometric alignment intensity values or a geometrically aligned data matrix. Figure 2 In This is the data matrix after alignment.
[0047] Wherein, time transformation function for:
[0048] .
[0049] in, Original three-dimensional feature tensor elements Align the time grid and preserve the time coordinates in the time axis. The number of time control points, For the first The position of each time control point on the timeline. , , Original three-dimensional feature tensor elements The transformation parameters of the learning time; This is a nonlinear kernel function. A nonlinear kernel function is a one-dimensional monotone radial basis function. One-dimensional monotone radial basis functions It must include at least one of the following forms:
[0050] 1. Cubic radial basis functions .
[0051] 2. Radial basis functions for thin plate splines ,in, ,Right now: To avoid positive decimals with unstable values.
[0052] 3. Linear combinations or numerically stable variants of cubic radial basis functions and thin-plate spline-type radial basis functions.
[0053] The time transformation function from the aligned time grid to the retained time axis is used to determine the time coordinates of the aligned reference time point in the retained time axis.
[0054] And, the time transformation parameters to be learned and It is a non-negative parameter; or, the parameter to be learned during time transformation. and Generation through activation functions to improve time transformation functions The global monotonicity.
[0055] Furthermore, when executing the nonlinear kernel function, the retained time axis is normalized to the interval [0,1] using linear or monotonic interpolation, and the position of the time control point on the retained time axis is set according to the following strategy:
[0056] A fixed set of time control points is formed by equidistantly placing points in the interval [0,1].
[0057] Alternatively, the position of the time control point in the [0,1] interval can be adaptively selected based on the total ion current peak density, gradient program change nodes, or prior knowledge.
[0058] The position of the time control point on the time axis is either a fixed parameter or can be fine-tuned within a limited range through programming with additional learnable offsets, in order to balance model expressiveness and training stability.
[0059] Time transformation function first derivative and second derivative Alignment time grid after discretization The above calculations are performed using finite differences or analytic derivatives. The first and second derivatives are used to construct the following regularization terms:
[0060] Smoothness regularization term:
[0061] .
[0062] Used to suppress high-frequency oscillations during time transformation and ensure smooth alignment curves.
[0063] Monotonicity regularization term:
[0064] .
[0065] This is used to penalize regions where the derivative is negative, thereby approximating the target value after training convergence. This ensures that the timeline is not reversed.
[0066] During the alignment process, the intensity transformation is defined as:
[0067] .
[0068] .
[0069] This is the formula for strength compensation. Wherein, Let be the intensity value of the i-th intermediate 3D feature tensor in the aligned time grid. To be based on the time transformation function The geometric alignment strength value obtained by interpolation sampling in the i-th original 3D feature tensor. Time transformation function The Jacobian factor, which is the time transformation function. The first derivative of the function is used as the intensity compensation factor. The intensity compensation factor is used to maintain the integral area (total amount of substance) of the chromatographic peak when a nonlinear scaling transformation occurs on the time axis, based on the principle of integral variable substitution. For a single chromatographic peak, its integral area remains essentially unchanged before and after the transformation. The principle of differentiable time transformation and peak area conservation compensation based on one-dimensional RBF is illustrated in the diagram below. Figure 4 As shown. Simultaneously, time transformation smoothness constraints and monotonicity constraints are introduced during training, and penalties are applied to regions where the first derivative is negative, in order to constrain... And it inhibits excessive distortion.
[0070] Peak area conservation is achieved through time sampling points The interpolated intensity is multiplied by the corresponding derivative. ,in Used to approximate peak area conservation, for a single chromatographic peak, its integral area remains basically unchanged before and after the transformation; in the discrete case, the peak area before and after the transformation is compared by numerical integration or trapezoidal formula to verify the effectiveness of the peak area conservation approximation, and a peak area difference penalty term can be added to the loss function as needed.
[0071] Among them, the micro-interpolation operators in interpolation sampling include one of linear interpolation, cubic spline interpolation, or piecewise polynomial interpolation, and interpolation is performed only along the retention time dimension.
[0072] Specifically, linear interpolation can be performed on adjacent sampling points along the retention time axis to ensure that the interpolation weights are continuously differentiable functions; or, under the premise of satisfying numerical stability, cubic spline interpolation or piecewise polynomial interpolation can be used; and the interpolation is only performed along the retention time dimension T, and not along the feature dimension F, thereby maintaining the semantic consistency and physical interpretability of the feature channels.
[0073] Step S130: Based on the sample batch identifier, injection order and acquisition time characteristics of the intermediate three-dimensional feature tensor, a scaling factor and a bias factor are generated through an explicit intensity drift correction layer. An element-wise affine transformation is then performed on the intermediate three-dimensional feature tensor based on the scaling factor and the bias factor to achieve intensity drift correction, thereby obtaining the target three-dimensional feature tensor. Figure 2 in, formula middle, This refers to the target 3D feature tensor output corresponding to the intensity drift correction layer. This is the scaling factor corresponding to the intensity drift correction layer; This is the bias factor corresponding to the intensity drift correction layer; The metadata corresponding to the intensity drift correction layer is the sample batch identifier (Batch ID) and the injection order.
[0074] This step achieves intensity drift correction for the intermediate 3D feature tensor through an explicit differentiable intensity drift correction layer. The input and output of the explicit differentiable intensity drift correction layer are the intermediate 3D feature tensor and the target 3D feature tensor, respectively. In the explicit differentiable intensity drift correction layer, firstly, based on the sample batch identifier, injection order, and acquisition time characteristics of the intermediate 3D feature tensor, a scaling factor and a bias factor are generated through a neural network mapping function. Secondly, based on the scaling factor and the bias factor, an element-wise affine transformation is performed on the intermediate 3D feature tensor to achieve intensity drift correction. The above principle is illustrated as follows: Figure 5 As shown.
[0075] The explicit intensity drift correction layer includes an embedded network, an intensity correction factor generation network, and an explicit correction layer.
[0076] Intensity drift correction is achieved by generating scaling and bias factors through an explicit intensity drift correction layer and performing element-wise affine transformations on the intermediate 3D feature tensor based on the scaling and bias factors. Specifically, this includes:
[0077] The embedding network is used to embed and map the sample batch identifier, injection order, and acquisition time of the intermediate three-dimensional feature tensor to obtain the corresponding vector representation.
[0078] The multiple vector representations of the intermediate three-dimensional feature tensor are concatenated or weighted and fused as input to the intensity correction factor generation network. The intensity correction factor generation network outputs the scaling factor and the bias factor. The scaling factor is activated by Softplus or exponential activation to ensure that the scaling factor is non-negative, and the bias factor is activated by linear or bounded activation to prevent numerical divergence. The intensity correction factor generation network includes a multilayer perceptron, which contains a nonlinear activation layer and a normalization layer to enhance the fitting ability to complex nonlinear drift patterns.
[0079] Intensity drift correction is achieved by performing an element-wise affine transformation on the intermediate three-dimensional feature tensor using the scaling factor and the bias factor through the explicit correction layer.
[0080] And, the formula for element-wise affine transformation is:
[0081] .
[0082] in, and They are respectively Scaling factor and bias factor, For the The transformed i-th target 3D feature tensor. Scaling and bias factors are physically valid through non-negative activation (e.g., Softplus) or parameter constraints. Since the scaling factor is a multiplicative factor, it should be constrained... (in It is a small positive number (such as 1e-6) to prevent division by zero errors or gradient vanishing.
[0083] Quality control stability and detrending constraints based on a reference anchor set are introduced into the loss function; the reference anchor set includes at least one of quality control samples, internal standards, or pseudo-reference samples based on in-batch biological sample statistics.
[0084] Step S140: The target three-dimensional feature tensor is enhanced by a dual-stream attention mechanism based on global flow features and differential flow features. The enhanced target three-dimensional feature tensor is then input into a quantitative regression module based on sequence modeling. The quantitative regression module outputs a quantitative feature matrix.
[0085] Feature enhancement of the target 3D feature tensor is achieved through a dual-flow attention mechanism based on global flow features and differential flow features. Specifically, this includes performing robust aggregation operations on the target 3D feature tensor along the sample batch dimension to obtain global flow features representing the commonality spectrum. And construct the difference flow (differential flow) features of each target 3D feature tensor relative to the global feature flow. A dual-stream attention mechanism is used to calculate attention weights through feature interaction, integrating global and differential stream features into the target 3D feature tensor to obtain enhanced target 3D features. The feature interaction calculation of attention weights utilizes G as the context for enhancement. That is, G acts as both a key and a value in self-attention. Used as a query value.
[0086] For example, the above process is as follows:
[0087] Robustly aggregate all target 3D feature tensors in the current micro-batch along the sample batch dimension to obtain the global reference spectrum. Robust aggregation includes at least one of the following: median, mean after removing extreme values, and weighted average.
[0088] For each target 3D feature tensor, through Construct differential flow features.
[0089] Differential flow characteristics Global streaming features are mapped to query vectors via linear projection. By mapping linear projection to keys and values, a scaled dot product attention mechanism is used to obtain the fused contextual features.
[0090] In the output stage, the context features and the target 3D feature tensor are fused through weighted summation, linear transformation after concatenation, or gating mechanism to enhance the separability of weakly different signals.
[0091] The enhanced target 3D feature tensor is then input into a quantitative regression module based on Transformer, Quanformer, or other sequence modeling networks. The quantitative regression module is configured to perform peak boundary refinement and response intensity normalization, and output a normalized quantitative feature matrix with rows as samples and columns as target substances or features. The output values include the corrected relative peak area, normalized intensity response, or (when external standards are available) the target substance concentration value. The output layer in the quantitative regression module uses a non-negative activation function to ensure that the predicted values are non-negative real numbers.
[0092] Specifically, the quantitative regression module is based on a Transformer-based feature refinement and integration architecture (such as the Quanformer model, from ZHANG Zhengyi, YANG Huan, WANG Yanyi, ZHANG Lei, LIN Shu-Hai. QuanFormer: A Transformer-Based Precise Peak Detection and Quantification Tool in LC-MS-Based Metabolomics[J]. Analytical Chemistry, 2025, 97(5):2698-2706. DOI:10.1021 / acs.analchem.4c04531.), which is configured to achieve quantification by predicting the integration boundary rather than directly regressing the concentration. The architecture specifically includes:
[0093] (1) Sequence embedding module: The time-feature two-dimensional structure is divided into patches according to the time axis or time-feature blocks, and each patch is mapped to an embedding vector through a linear layer, and position encoding is superimposed to form an input sequence containing temporal information.
[0094] (2) Global Dependency Modeling Module: Stack at least two layers of multi-head self-attention and feedforward network modules to perform global context encoding on the input sequence to capture long-range baseline drift and chromatographic peak distortion features.
[0095] (3) Boundary Detection and Mask Generation Module: At the end of the network, instead of directly outputting a single quantitative value, the decoding head predicts the start and end integral boundary probability distribution of the target substance on the retention time axis or the binary mask of the region of interest (ROI). The principle of boundary detection and integral quantification based on hybrid supervision is illustrated in the diagram below. Figure 6 As shown.
[0096] (4) Differentiable Integral Quantitative Module: Based on the predicted boundary or mask, the response intensity of the target substance is calculated on the characteristic tensor after intensity drift correction by differentiable weighted summation or integral operator, thereby outputting a standardized quantitative value with physical interpretability.
[0097] For the joint loss function, the joint loss function simultaneously optimizes all or part of the parameters of the one-dimensional monotonic radial basis function space transformation network, explicit intensity drift correction layer, two-stream attention mechanism, and quantitative regression module in the same computational graph through a single backpropagation. Through a virtualized columnar data access interface, the final alignment parameters, intensity correction factors, and quantitative results are incrementally written to the storage medium in a columnar manner, and the global logical index is updated to support incremental inference and fine-tuning of subsequent batches of data.
[0098] Joint loss function for:
[0099] .
[0100] .
[0101] .
[0102] .
[0103] .
[0104] .
[0105] in, , , , , and These are, respectively, the quantitative regression error term, the retention time alignment consistency constraint, the quality control stability and detrending constraint, the time transformation monotonicity regularization term, the time transformation smoothness regularization term, and the L2 norm regularization constraint on the network parameters. , , , , and All are non-negative weighting coefficients. For a pre-defined set of quality control (QC) target substances or characteristics, It is a dispersion measure operator used to calculate variance, standard deviation, or coefficient of variation (CV). This is a correlation measurement operator used to calculate the Pearson correlation coefficient, Spearman correlation coefficient, or mutual information. Let f be the intensity value of the i-th sample in the complete time series after intensity correction on the feature channel f. This is the time-averaged intensity or peak area response value of sample i on the characteristic channel f after correction. For QC sample set; The degree of dispersion; To determine the correlation with injection order or collection time.
[0106] The quantitative regression error term includes label-based supervised loss and self-supervised consistency loss.
[0107] Quality control stability and detrending constraints include one or a combination of the following: statistical dispersion minimization constraint, temporal correlation decoupling constraint, and adaptive baseline constraint.
[0108] The quantitative regression error term specifically includes: (a) Reference Feature Supervision Loss: supervising the model output using predefined reference quantitative targets; the reference quantitative targets include regions of interest (ROI) boundaries, peak area integrals, or pseudo-labels derived from expert annotations or high-confidence results extracted by traditional algorithms; (b) Self-supervised Consistency Loss: minimizing the coefficient of variation (CV) or manifold distribution differences of the same substance among different technical replicas or quality control samples in the absence of external concentration labels; and (c) (optional) absolute concentration regression error when standards are available.
[0109] Quality control sample stability / detrend constraint Configured as a data harmonization and representation alignment objective function to ensure that the feature representations learned by the network maintain statistical distribution consistency throughout the entire collection period of quality control samples, thereby suppressing non-biological batch effects; constraints include at least one of the following or a combination thereof:
[0110] (1) Statistical dispersion minimization constraint: For a preset target substance or feature set, calculate at least one of the following dispersion measures of the strength of the quality control sample after correction: variance, standard deviation, coefficient of variation (CV), and mean absolute deviation (MAD), and minimize the dispersion measure to achieve response stability across time and / or batch.
[0111] (2) Temporal correlation decoupling constraint: Calculate at least one of the following correlation measures, namely Pearson correlation coefficient, Spearman correlation coefficient or mutual information, for the change of quality control sample strength with injection order, collection time or operation number, and minimize the absolute value, square or weighted term of the correlation measure to weaken the time-related system drift trend.
[0112] (3) Alternative strategy: Adaptive baseline constraint: In the absence of external quality control samples or when the quality control samples are significantly different from the biological sample matrix, making quality control constraints inapplicable, a robust statistical baseline based on the biological sample population (including at least one of the sliding window median, quantile baseline or low-rank approximation background) is constructed, and the constraint objective is configured to minimize the systematic shift of the biological samples relative to the statistical baseline, as an alternative or supplement to the quality control stability / detrending constraint.
[0113] Weight coefficients of each term in the joint loss function The setting methods include: adjusting within a preset range based on experience so that the various losses are roughly equal in magnitude in the early stage of training; or using uncertainty weighting, adaptive balancing, or gradient norm-based balancing methods to dynamically update the weights of various items during training, so as to avoid a certain loss item dominating the optimization process or some regularization items being ignored.
[0114] When performing inference or incremental training on new batches of data, the following strategies are adopted: use global logical indexes to locate the trained model and its corresponding alignment and intensity correction parameter versions; in inference-only mode, freeze all network parameters and perform only forward computation to obtain the alignment results and quantitative matrix of the new samples; in incremental fine-tuning mode, freeze some backbone network layers related to peak shape modeling and fine-tune only a few parameters in the time transformation layer and the subnetwork related to intensity drift correction to quickly adapt to the new batch effect while maintaining the comparability and traceability of historical results.
[0115] As explained above, this embodiment employs an end-to-end deep learning network to analyze data from a large-scale chromatography-mass spectrometry cohort. During the unified backpropagation process of the end-to-end deep learning network, retention time alignment, cross-batch intensity drift correction, and quantitative regression are performed simultaneously. This means that retention time alignment, cross-batch intensity drift correction, and quantitative regression are considered as a unified joint operation, rather than sequential operations, overcoming the shortcomings of the step-by-step strategies in traditional techniques. Specifically, retention time alignment, cross-batch intensity drift correction, and quantitative regression are implemented through different network modules within the end-to-end deep learning network. Retention time alignment is achieved through a one-dimensional monotonic radial basis function space transformation network; cross-batch intensity drift correction is achieved through an explicit differentiable intensity drift correction layer; and quantitative regression is achieved through a two-stream attention mechanism and a quantitative regression module. The network parameters of these modules are all synchronously optimized based on a joint loss function.
[0116] Specifically, this embodiment has the following technical effects:
[0117] 1. I / O Bottleneck Breakthrough and Universality (IO Virtualization): This embodiment innovatively adopts a virtualized columnar data access interface, supporting direct indexing and logical columnar slicing of raw heterogeneous mass spectrometry data (such as .raw, .mzML, etc.) without mandatory pre-conversion to a specific format. Furthermore, when data is pre-converted to a columnar or block storage format, higher throughput can be achieved by utilizing built-in indexes and block compression. Combining streaming micro-batch loading and memory buffering technology, it effectively solves the risk of memory overflow (OOM) and I / O throughput bottlenecks in large-scale queue (thousands of cases or more) analysis, significantly reducing the time cost and storage space occupation of data preprocessing.
[0118] 2. End-to-End Error Cooperative Suppression: This approach abandons the traditional step-by-step error accumulation model of "alignment first, then correction, then quantification," and constructs an end-to-end integrated deep learning framework that includes retention time alignment, intensity drift correction, and quantitative regression. By constructing a joint loss function, the alignment and correction parameters at the front end are directly supervised and fed back by the quantitative accuracy at the back end, achieving global cooperative optimization of all parameters across the entire process.
[0119] 3. Dual Adaptive Calibration: To address complex instrument fluctuations in large-scale acquisition queues, this embodiment incorporates a differentiable alignment layer of a one-dimensional monotonic radial basis function (RBF) spatial transformation network (1D-RBF-STN) and an explicit intensity drift correction layer. The former adaptively calibrates the nonlinear warping of retention time, while the latter accurately simulates and removes the nonlinear signal decay of the instrument based on injection sequence and batch information, thereby ensuring extremely high quantitative consistency in complex matrices and long-cycle acquisition.
[0120] 4. High-sensitivity signal detection (Signal Enhancement): Introducing a unique global-differential dual-stream attention mechanism, which uses a globally stable common background as a reference context, and effectively amplifies weak sample-specific difference signals through feature interaction, significantly improving the detection rate and signal-to-noise ratio of low-abundance disease biomarkers under complex background noise.
[0121] 5. Industrial-grade high-throughput parallelism and standardized output: This embodiment employs a naturally compatible parallel computing architecture with multi-core CPUs and multi-GPUs. Micro-batch quantization processing enables efficient utilization of computing resources, significantly improving analysis speed. Simultaneously, the system directly outputs a standardized quantitative feature matrix with samples as rows and target substances as columns, eliminating the cumbersome intermediate format conversion steps and seamlessly connecting to downstream statistical analysis and biomarker discovery workflows.
[0122] The following is an example workflow for a joint alignment and quantification method for large-scale mass spectrometry queues.
[0123] This example runs on a high-performance computing platform, such as a multi-socket server CPU with 64 cores or more, equipped with 4–8 GPUs (≥24GB of video memory), and the system software environment can be a Linux operating system, and the deep learning framework can be PyTorch, TensorFlow or a similar platform.
[0124] The overall process includes the following steps:
[0125] S1. Construction of virtualized columnar interface and micro-batch quantization.
[0126] (1) Metadata parsing and index construction.
[0127] The system receives the input set of raw chromatography-mass spectrometry data files. File formats can include, but are not limited to, vendor-specific formats (such as .raw, .wiff, .d) and standard open formats (such as .mzML). To avoid memory and I / O bottlenecks caused by decoding all spectral data at once, this example first performs a lightweight metadata parsing step, including:
[0128] Read the file header, index list, offset table, or directory tree structure; extract the retention time array, scan number, total ion current (TIC) summary, and corresponding physical file pointer (byte start offset, byte length, etc.) for each scan; organize the above metadata into a memory-resident index structure, such as a hash table or B+ tree, where the key can be (sample ID, retention time window ID) and the value can be (file path, start offset, data block size, mass-to-charge ratio range information, etc.).
[0129] At this stage, the complete m / z-intensity matrix is not read; only the metadata is processed, resulting in fast parsing speed and low memory usage.
[0130] (2) Virtualized columnar access interface.
[0131] Based on the aforementioned metadata index, this example constructs a virtualized columnar data access interface. This interface exposes a logical view of "feature / column access" to the upper-layer deep learning modules, while the underlying actual I / O operations still target the original row-oriented storage files. The virtualization interface is responsible for: calculating the starting offset and length of the corresponding data block on disk based on the user-requested sample index, retention time interval, m / z interest range, etc.; supporting parallel fseek / fread operations to read micro-batch segments from multiple files as needed; and for data files pre-converted to columnar or block formats such as HDF5, Zarr, Lance, and Parquet, the virtualization interface can directly call the columnar read APIs for these formats, combining block compression and indexing to improve throughput.
[0132] This virtualized columnar access interface decouples the organization of the original physical files from the upper-layer computing logic, enabling unified access to heterogeneous data sources.
[0133] (3) Streaming micro-batch loading and three-dimensional tensor quantization.
[0134] During the training or inference phase, this embodiment uses a producer-consumer model for streaming data loading: a prefetch thread pool is maintained on the CPU side as the "producer"; the deep learning training process acts as the "consumer", retrieving micro-batch data from the queue on demand.
[0135] When forming each microbatch: select (B) samples from the global sample list according to the predetermined microbatch size (B); for each sample, binning or feature extraction is performed in the m / z dimension within a specified retention time range (e.g., truncated or resampled to a uniform time grid of length (T)) to construct a three-dimensional feature tensor. ,in, For sample batch dimension, To retain the time dimension for discrete elements, This refers to the mass-to-charge ratio or the dimension of the pre-extracted features.
[0136] To reduce memory pressure, this embodiment employs a sparse coding strategy: during the data loading stage, a threshold is set for the signal strength, retaining only entries that are non-zero or higher than the noise level, and recording ((b,t,f)) and their corresponding strength values in coordinate (COO) format. Subsequently, before entering the network segment that requires convolution or attention operations, the sparse format is converted into locally dense blocks or moderately compressed dense tensors to balance memory utilization and operator compatibility.
[0137] Meanwhile, the system can enable a memory-aware scheduling mechanism: monitor GPU memory usage in real time, and when it approaches a preset threshold, automatically reduce the micro-batch size, shorten the time window, or sparsify the standard to avoid OOM (Out-Of-Memory) errors.
[0138] S2. Alignment of time transformation and differentiable interpolation of one-dimensional monotonic radial basis functions.
[0139] In a preferred embodiment, the present invention models retention time alignment as a trainable and deployable "differentiable time transformation layer" in an end-to-end network, which is collaboratively completed by a localization sub-network, a grid generator, and a differentiable sampler. The localization sub-network extracts drift-related temporal features (e.g., the total ion flow curve and its peak distribution characteristics obtained from the aggregation of three-dimensional feature tensors) from the observational information of the sample to be aligned, and outputs a set of parameters characterizing the time transformation function of that sample, allowing different samples to obtain their own adaptive time transformations. The grid generator receives the aforementioned time transformation parameters and a unified reference time grid, generating mapping coordinates from the reference time axis to the original time axis, thereby determining which retention time position in the original data each alignment time point should be valued at. The differentiable sampler performs interpolation resampling in the retention time dimension based on the mapping coordinates, outputting a geometrically aligned feature or intensity matrix; in a preferred embodiment, the resampling intensity can also be compensated by combining the local scaling of the time mapping to reduce the impact of time stretching / compression on peak area and quantitative consistency.
[0140] During the training phase, this invention integrates alignment consistency constraints, quality control stability constraints, and quantitative regression errors into a joint loss function. Utilizing a chained gradient backpropagation mechanism based on differentiable sampling and differentiable mapping, error signals generated by downstream quantitative and quality control targets can be propagated back to the time mapping parameters, further updating the relevant parameters of the localization sub-network and alignment module. This achieves collaborative optimization of "alignment-correction-quantification" within the same computational graph. During the inference phase, the localization sub-network can directly provide the time transformation parameters and complete the alignment for a new sample in a single forward computation, eliminating the need for iteratively solving the time transformation parameters for each new sample. Compared to directly treating the time transformation parameters as free variables for each sample and repeatedly optimizing them during inference, this invention, through its "parameter prediction amortization inference + differentiable resampling" structural design, significantly reduces computational and storage overhead in large-scale queue scenarios, improving the ability to quickly access new samples and generalize. Simultaneously, through structural constraints such as monotonicity and smoothness, and interpolation resampling mechanisms, it reduces the risk of unreasonable distortions in the alignment process, improving training stability and cross-batch quantitative consistency. This, in turn, reduces the accumulation of errors in multi-stage serial processes and enhances the reliability of overall qualitative and quantitative output.
[0141] To address the nonlinear drift of chromatographic retention times between different samples, this example introduces a one-dimensional monotonic radial basis function (RBF) time transformation layer at the front end of the neural network and designs it as a differentiable structure to support end-to-end training.
[0142] Time transformation function definition: For the i-th sample in the microbatch, this embodiment defines a time transformation function on the normalized retention time interval ([0,1]):
[0143] .
[0144] in: ; Set the preset number of time control points, for example, 8–32; For the first The location of each control point can be determined by equidistant points or by prior experience. For the sample The parameters to be learned; For radial basis functions, it can be: , Or its numerically stable variant. The above form is similar to a one-dimensional radial basis function interpolation model, combining flexibility and controllability.
[0145] Localization Subnetwork and Parameter Prediction: To make the time transformation function adaptable to each sample, this example constructs a localization subnetwork as the time transformation layer to predict the parameter set.
[0146] .
[0147] The input to the localization subnetwork is the sample. The spectral features can be derived from the three-dimensional feature tensor. Obtained through the following steps:
[0148] Along the feature dimension Aggregation (summation, maximum value, or weighted average) yields the total ion current (TIC) time series. .
[0149] Will Input several layers of one-dimensional convolutional networks (1D Conv) or lightweight Transformers to extract peak shapes and temporal patterns.
[0150] The aforementioned intermediate features are then subjected to global pooling (average pooling or max pooling) along the time dimension to obtain a fixed-length vector. .
[0151] Subsequently, a multilayer perceptron (MLP) was used to... Mapped to parameter vector Furthermore, activation functions such as Softplus and Sigmoid can be used to impose non-negative or bounded constraints on some parameters, which is beneficial for monotonicity and numerical stability.
[0152] From the above form, we can obtain the first derivative:
[0153] .
[0154] And the second derivative:
[0155] .
[0156] In actual implementation, (r), (r) can be analytically calculated based on the selected radial basis function and in a discrete-time grid. Numerical approximation is performed on the above.
[0157] Differentiable interpolation sampling and time alignment:
[0158] Obtaining the time transformation function This example then uses a unified reference time grid. Upsampled samples The three-dimensional feature tensor. Let the features before transformation be denoted as . The alignment features after transformation are Defined as:
[0159] .
[0160] in, It is obtained by performing linear or higher-order interpolation on adjacent discrete points on the original time axis; A local scale characterizing time stretching is used to compensate for peak area during time compression or stretching, in order to approximately maintain peak area conservation.
[0161] Due to interpolation weights and Both the function and its derivative are differentiable functions, and the overall sampling process is differentiable with respect to the upstream parameters of the network, allowing for backpropagation through an automatic differentiation framework.
[0162] Smoothness and monotonicity regularization constraints:
[0163] To prevent unreasonable distortion of the timeline and maintain the timing logic of RT, this example introduces the following regularization term during training:
[0164] Smoothness constraint:
[0165] .
[0166] Used to suppress high-frequency oscillations during time transformation and ensure smooth alignment curves.
[0167] Monotonicity regularization term:
[0168] .
[0169] Penalize the region where the first derivative is negative to approximate the region. The condition for monotonically increasing.
[0170] The regularization term mentioned above will be optimized together with other terms in the subsequent joint loss function through weighted combination.
[0171] S3. Explicit Differentiable Strength Drift Correction and QC Constraints.
[0172] After completing the retention time alignment, this embodiment further performs explicit modeling and correction of intensity drift across batches.
[0173] Correction factor generation network:
[0174] For sample i, construct its metadata feature vector:
[0175] .
[0176] in: Indicates the batch number to which the sample belongs; This indicates the injection order, which can be normalized to ([0,1]); It may include features such as the date and time of collection and the operation number; The embedding function is for discrete identifiers. Metadata features are input into a multilayer perceptron (MLP):
[0177] .
[0178] And further map to generate intensity correction factors:
[0179] .
[0180] .
[0181] in and It can be related to the feature dimension A compatible one-dimensional vector is applied to all time steps via a broadcast mechanism; scaling factors are guaranteed using Softplus or similar non-negative activation. . and All are mapping coefficients; , All are mapping constants.
[0182] Intensity Affine Correction:
[0183] Features that are already aligned to RT Perform an explicit affine transformation:
[0184] .
[0185] parameter The main description focuses on the multiplicative decay of instrument sensitivity over time or batches. This simulates additive offsets such as baseline drift. Through end-to-end training, the network can automatically learn systematic intensity variation patterns across batches.
[0186] QC Sample Stability and Detrending Constraints:
[0187] In large-scale cohorts, quality control (QC) samples are typically interspersed. This embodiment utilizes QC samples to apply supervised constraints to the intensity correction effect. Let the QC sample set be... For a specific characteristic channel or target substance, the corrected intensity sequence on the QC sample is calculated as follows:
[0188] Dispersion (e.g., variance, standard deviation, coefficient of variation):
[0189] .
[0190] Correlation with injection order or acquisition time :
[0191] .
[0192] in Can be a sample In the passage The time-averaged intensity or peak area.
[0193] By minimizing and This can force the intensity correction layer to remove systematic drift that is highly correlated with time / sequence, thereby improving the stability and comparability of the corrected data.
[0194] Adaptive baseline constraints (no QC scenario):
[0195] For scenarios where physical QC samples may be lacking or the matrix effect of QC samples may not be applicable in large clinical cohorts, this example also configures an adaptive baseline mode as a backup strategy.
[0196] In this mode, the system no longer relies on predefined... Instead of using a collection, it utilizes the collection of all biological samples within the current microbatch or sliding window. Construct a virtual reference baseline. The specific steps are as follows:
[0197] 1. Trendline Construction: For each feature channel Calculate the corrected intensity of biological samples Local trend terms with injection sequence or time This trend can be measured by... The systematic background drift at the population level is obtained by performing a sliding median filter with a large window width or by low-order polynomial fitting.
[0198] 2. Optimize Objective Restructuring: Shift the constraint objective from "minimizing the dispersion of individual variables" to "minimizing the fluctuation of the trend term." Construct a loss function to penalize... The variance or its correlation with the collection time .
[0199] By employing the above strategy, the model eliminates only the systematic drift (detrending) at the population level, while allowing individual samples to retain their inherent discrete distribution above and below the trend line. This achieves batch effect correction without the need for external standards, while avoiding the loss of biological differences due to overcorrection.
[0200] S4. Global-Differential Two-Stream Attention Enhancement Module.
[0201] To highlight the true differences between samples amidst complex background noise, this embodiment designs a global-differential dual-stream attention module to further extract and enhance the intensity-corrected features.
[0202] Common flow feature construction:
[0203] For a micro-batch, let the intensity-corrected characteristic be... Along the sample dimension Perform robust aggregation to obtain global common features. :
[0204] .
[0205] in, It can be the mean, median, or weighted average after removing extreme values, used to represent the reference peak shape of the material in the current batch on the aligned time axis.
[0206] Differential flow feature construction:
[0207] For each sample Define its residual or difference characteristics:
[0208] .
[0209] This residual feature characterizes the degree of deviation of a sample from the common background spectrum, and is more likely to reflect the biological differences caused by disease state or treatment conditions.
[0210] Dual-stream attention fusion:
[0211] Will and Treating it as a two-dimensional feature map, it can first be flattened into a sequence along the time-feature dimension, and then the attention input can be constructed:
[0212] Differential flow characteristics Perform a linear transformation to obtain the query vector :
[0213] .
[0214] ( ) represents a linear transformation.
[0215] Obtain the key by performing a linear transformation on the common flow characteristics. Sum :
[0216] .
[0217] .
[0218] , as well as Both represent linear transformation weights; , as well as All represent constants corresponding to linear transformations, and G represents the global flow characteristic.
[0219] By scaling the dot product attention:
[0220] .
[0221] in This is the attention feature dimension. Then... Reshape back to the time-feature shape and compare with the original difference features. By splicing or weighted fusion, enhanced features are obtained:
[0222] .
[0223] in, It can be a weighted sum, a concatenated result after linear layer mapping, or a gating mechanism.
[0224] Through the aforementioned dual-stream attention structure, this example utilizes the common stream to provide a stable background context, making it easier for the network to capture weak signals, low-abundance features, and local deformations on the differential stream, significantly improving its sensitivity to subtle biological effects.
[0225] S5. Sequence Quantitative Regression Network and Non-negative Output.
[0226] Features after attention enhancement As input to a quantitative regression network. This example uses a Transformer-based sequence regression structure such as Quanformer, but CNN, RNN, or a combination thereof can also be used.
[0227] Serialization and Embedding:
[0228] Will Patch segmentation is performed in the time dimension or the time-feature two-dimensional dimension, for example, every certain number of time steps is a patch, and then it is mapped to a hidden vector sequence through a linear layer or low-dimensional convolution, and then the positional encoding is superimposed and used as the input of the Transformer encoder.
[0229] Transformer encoder:
[0230] Several layers of self-attention-feedforward network modules are stacked to capture long-range dependencies between different time periods or peaks through multi-head self-attention, and training stability is maintained through residual connections and layer normalization.
[0231] Quantitative regression output layer:
[0232] The output of the Transformer is pooled (e.g., using the [CLS] token or global average pooling), and then mapped to a dimension of [dimensionality missing] through several fully connected layers. The vector, where The target substance or feature quantity. Finally, a non-negative activation function such as Softplus, ReLU, or Sigmoid-Scale is used.
[0233] .
[0234] To ensure the predicted quantitative value It is a non-negative real number.
[0235] S6. Joint Loss Function and End-to-End Optimization.
[0236] To achieve integrated optimization of alignment, intensity correction, and quantitative regression, this example constructs the following joint loss function:
[0237] .
[0238] in:
[0239] : Quantitative regression error term, which can be the mean squared error (MSE), Huber loss, or mean absolute error (MAE).
[0240] : Retain time alignment consistency constraints, which can be constructed based on the similarity of peak shapes after alignment, the correlation coefficient between samples, etc.
[0241] QC stability and detrend constraint, including the aforementioned , .
[0242] Time transformation monotonicity regularization; Time transformation smoothness regularization; Regularization constraints such as the L2 norm of network parameters; These are non-negative weighting coefficients, which can be adjusted through empirical settings or adaptive strategies.
[0243] By calculating the loss on the same computation graph and performing a backpropagation, all trainable parameters in the localization subnetwork, time transformation layer, intensity correction layer, dual-stream attention module, and quantitative regression network can be updated synchronously, achieving end-to-end collaborative optimization and avoiding the error accumulation problem of traditional multi-step processes.
[0244] S7. Global logical index and incremental update mechanism.
[0245] To support continuous analysis of large-scale queues and incremental access to new batches of data, this example maintains a global logical index structure to record the metadata of each sample and the location of its analysis results. This index includes at least:
[0246] Sample identifier, batch identifier, injection order; currently used alignment parameter version number, intensity correction factor version number; pointer to the position of the quantitative result in columnar storage, etc.
[0247] When new batches of data are added, the system can perform forward inference or minor incremental fine-tuning only on the new samples. In one implementation, the backbone network layers related to peak shape modeling can be frozen, and only the subnetwork parameters related to intensity drift correction can be updated to ensure the traceability of historical results and the consistency of quantification across batches, while avoiding duplicate calculations of all historical data.
[0248] S8. Model training strategy.
[0249] The reference quantitative target in the quantitative regression error term can be flexibly selected based on data acquisition costs: In a preferred implementation (high-performance mode), a strategy similar to Quanformer is adopted, using peak boundaries and integral areas verified by human experts as strong supervision signals to train the network to learn fine chromatographic peak morphology features; In another implementation (high-throughput / low-cost mode), a weak supervision strategy is adopted, directly using peak areas extracted from high signal-to-noise ratio samples by traditional algorithms (such as XCMS, OpenMS) or software as pseudo-labels, i.e., rough boxes, to drive the network to approximate the integral performance of traditional algorithms through large-scale data and to suppress integral noise of individual samples by utilizing the smoothing characteristics of deep learning.
[0250] In one exemplary embodiment, a computer device is also provided, including a memory and a processor, wherein the memory stores a computer program, and the processor executes the computer program to implement the steps in the above-described method embodiments.
[0251] In one exemplary embodiment, a computer-readable storage medium is provided storing a computer program that, when executed by a processor, implements the steps in the above-described method embodiments.
[0252] In one exemplary embodiment, a computer program product is provided, including a computer program that, when executed by a processor, implements the steps in the above-described method embodiments.
[0253] It should be noted that the user information (including but not limited to user device information, user personal information, etc.) and data (including but not limited to data used for analysis, data stored, data displayed, etc.) involved in this application are all information and data authorized by the user or fully authorized by all parties, and the collection, use and processing of the relevant data must comply with relevant regulations.
[0254] Those skilled in the art will understand that all or part of the processes in the above embodiments can be implemented by a computer program instructing related hardware. The computer program can be stored in a non-volatile computer-readable storage medium. When executed, the computer program can include the processes of the embodiments described above. Any references to memory, databases, or other media used in the embodiments provided in this application can include at least one of non-volatile and volatile memory. Non-volatile memory can include read-only memory (ROM), magnetic tape, floppy disk, flash memory, optical memory, high-density embedded non-volatile memory, resistive random access memory (ReRAM), magnetic random access memory (MRAM), ferroelectric random access memory (FRAM), phase change memory (PCM), graphene memory, etc. Volatile memory can include random access memory (RAM) or external cache memory, etc. By way of illustration and not limitation, RAM can take many forms, such as Static Random Access Memory (SRAM) or Dynamic Random Access Memory (DRAM).
[0255] The databases involved in the embodiments provided in this application may include at least one type of relational database and non-relational database. Non-relational databases may include, but are not limited to, blockchain-based distributed databases. The processors involved in the embodiments provided in this application may be general-purpose processors, central processing units, graphics processing units, digital signal processors, programmable logic devices, quantum computing-based data processing logic devices, etc., and are not limited to these.
[0256] The technical features of the above embodiments can be combined in any way. For the sake of brevity, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, they should be considered to be within the scope of this specification.
[0257] This document uses specific examples to illustrate the principles and implementation methods of this application. The descriptions of the above embodiments are only for the purpose of helping to understand the methods and core ideas of this application. Furthermore, those skilled in the art will recognize that, based on the ideas of this application, there will be changes in the specific implementation methods and application scope. Therefore, the content of this specification should not be construed as a limitation of this application.
Claims
1. A joint alignment and quantification method for large-scale mass spectrometry cohorts, characterized in that, A quantitative feature matrix is output through an end-to-end deep learning network; wherein, during the unified backpropagation process of the end-to-end deep learning network, retention time alignment, cross-batch intensity drift correction, and quantitative regression are simultaneously optimized based on a joint loss function. The training steps of the end-to-end deep learning network include: The original three-dimensional feature tensors of all batches are obtained. The original three-dimensional feature tensors of one batch are used as standard samples, and the original three-dimensional feature tensors of the remaining batches are aligned with the standard samples. The original three-dimensional feature tensors are extracted from the chromatography-mass spectrometry queue. A one-dimensional monotonic radial basis function space transformation network is used to learn a time transformation function for each sample, from aligning the time grid to preserving the time axis, to achieve time alignment of each original three-dimensional feature tensor. Intensity compensation is performed on the original three-dimensional feature tensor based on the time transformation function and its first derivative. The intensity compensation is used to maintain the conservation of the chromatographic peak integral area when the time axis is nonlinearly stretched. The original three-dimensional feature tensor after time alignment and intensity compensation is the intermediate three-dimensional feature tensor. Based on the sample batch identifier, injection order, and acquisition time of the intermediate three-dimensional feature tensor, a scaling factor and a bias factor are generated through an explicit intensity drift correction layer, and an element-wise affine transformation is performed on the intermediate three-dimensional feature tensor based on the scaling factor and the bias factor to achieve intensity drift correction, thereby obtaining the target three-dimensional feature tensor. The target 3D feature tensor is enhanced by a dual-stream attention mechanism based on global flow features and differential flow features. The enhanced target 3D feature tensor is then input into a quantitative regression module based on sequence modeling, which outputs the quantitative feature matrix.
2. The joint alignment and quantification method for large-scale mass spectrometry cohorts according to claim 1, characterized in that, Obtain the original 3D feature tensors for all batches, specifically including: Parse the metadata index of the raw mass spectrometry data file, and establish an index mapping relationship from physical storage address to logical feature space based on the metadata index; Without decoding all the original mass spectrometry data files, data blocks in heterogeneous row or column storage files on the disk are read on demand through the index mapping relationship, and multiple data blocks are assembled into a batch of the original three-dimensional feature tensors in a streaming manner. The original three-dimensional feature tensor has dimensions including sample batch dimension, discrete retention time dimension, and feature dimension; The definition of the feature dimension includes: The preset mass-to-charge ratio range is divided into multiple non-overlapping or partially overlapping sub-intervals, and each sub-interval corresponds to one or more indices in the feature dimension; Alternatively, the pre-extracted spectral features can be mapped to a set of channels along the feature dimension, wherein the pre-extracted spectral features include pre-identified regions of interest, specific isotope distribution clusters, and fragment ion clusters; Elements in constructing the original three-dimensional feature tensor At that time, the data block is placed in the time window And fall into the feature unit f One or more ion signal intensities within the corresponding mass-to-charge ratio range are aggregated to obtain aggregated intensity values. These aggregated intensity values are then subjected to logarithmic transformation, range normalization, or standard deviation standardization to improve the comparability of different original three-dimensional feature tensors along the feature dimension and suppress the impact of extreme intensity values on training. The processed aggregated intensity values are then used as the element. During the aggregation process, the aggregation operators used include at least one of summation, maximum value taking, weighted average, and peak area integration.
3. The joint alignment and quantification method for large-scale mass spectrometry cohorts according to claim 1, characterized in that, The time transformation function for: ); in, Original three-dimensional feature tensor elements Align the time grid and preserve the time coordinates in the time axis. The number of time control points, For the first The position of each time control point on the timeline. , , Original three-dimensional feature tensor elements The transformation parameters of the learning time; It is a nonlinear kernel function; The one-dimensional monotonic radial basis function space transformation network includes a localization subnetwork, a mesh generator, and a differentiable sampler. The localization subnetwork, based on the input original three-dimensional feature tensor, outputs transformation parameters characterizing the time transformation function. These transformation parameters include a set of time control points defined in the normalized time domain [0,1]. The positions of these time control points are adaptively determined based on the peak density and / or gradient program variation nodes and / or prior anchor points of the total ion current curve. These time control points are either fixed parameters or finely adjusted within a preset range using learnable offsets while maintaining monotonicity. The mesh generator generates a mapping coordinate between an aligned time mesh and the original retained time axis based on the transformation parameters. The differentiable sampler performs differentiable interpolation sampling on the original three-dimensional feature tensor in the retained time dimension based on the mapping coordinates to obtain geometric alignment strength values or a geometrically aligned data matrix. When executing the nonlinear kernel function, the retained time axis is normalized to the interval [0,1] using linear or monotonic interpolation, and the position of the time control point on the retained time axis is set according to the following strategy: A fixed set of time control points is formed by equidistantly placing points in the [0,1] interval; Alternatively, the position of the time control point in the [0,1] interval can be adaptively selected based on the total ion current peak density, gradient program change nodes, or prior knowledge. The position of the time control point on the retained time axis is a fixed parameter or can be fine-tuned within a limited range through programming with additional learnable offsets, so as to balance the model's expressive power and training stability. The nonlinear kernel function is a one-dimensional monotonic radial basis function. The one-dimensional monotonic radial basis function It must include at least one of the following forms: cubic radial basis functions ; The Euclidean distance between the time control points; Thin plate spline type radial basis function ,in, To avoid positive decimals with unstable values; A linear combination or numerically stable variant of the cubic radial basis function and the thin-plate spline type radial basis function; And, the time transformation parameters to be learned and It is a non-negative parameter; or, the time transformation parameter to be learned and Generation through activation functions to improve time transformation functions The global monotonicity; The time transformation function first derivative and second derivative After discretization The above calculations are performed using finite differences or analytic derivatives.
4. The joint alignment and quantification method for large-scale mass spectrometry cohorts according to claim 3, characterized in that, The formula for the strength compensation is: in, Let i be the intensity value of the i-th intermediate 3D feature tensor in the aligned time grid. To be based on the time transformation function The geometric alignment strength value obtained by interpolation sampling in the i-th original three-dimensional feature tensor. Time transformation function The Jacobian factor, which is the time transformation function. The first derivative; The micro-interpolation operators in the interpolation sampling include one of linear interpolation, cubic spline interpolation, or piecewise polynomial interpolation, and the interpolation is performed only along the retention time dimension.
5. The joint alignment and quantification method for large-scale mass spectrometry cohorts according to claim 1, characterized in that, The explicit intensity drift correction layer includes an embedded network, an intensity correction factor generation network, and an explicit correction layer; Intensity drift correction is achieved by generating scaling and bias factors through an explicit intensity drift correction layer and performing element-wise affine transformations on the intermediate 3D feature tensor based on the scaling and bias factors. Specifically, this includes: The intermediate three-dimensional feature tensor is embedded and mapped using the embedding network to obtain the corresponding vector representation; The multiple vector representations of the intermediate 3D feature tensor are concatenated or weighted and fused as input to the intensity correction factor generation network. The intensity correction factor generation network outputs the scaling factor and the bias factor. The scaling factor is activated by Softplus or exponential activation to ensure that the scaling factor is non-negative, and the bias factor is activated by linear or bounded activation to prevent numerical divergence. The intensity correction factor generation network includes a multilayer perceptron, which contains a nonlinear activation layer and a normalization layer to enhance the fitting ability to complex nonlinear drift patterns. Intensity drift correction is achieved by performing an element-wise affine transformation on the intermediate three-dimensional feature tensor using the scaling factor and the bias factor through the explicit correction layer. The formula for the element-wise affine transformation is: ; in, and They are respectively Scaling factor and bias factor, For the The transformed i-th target three-dimensional feature tensor.
6. The joint alignment and quantification method for large-scale mass spectrometry cohorts according to claim 1, characterized in that, Feature enhancement of the target 3D feature tensor is performed using a dual-stream attention mechanism based on global and differential flow features, specifically including: A robust aggregation operation is performed on the target 3D feature tensor along the sample batch dimension to obtain global flow features representing the commonality spectrum. And construct the three-dimensional feature tensor of each target. Differential flow features relative to global feature flow ; The dual-stream attention mechanism calculates attention weights through feature interaction to integrate the global stream features and the differential stream features into the target 3D feature tensor, resulting in enhanced target 3D features. The feature interaction calculation of attention weights utilizes G as context for enhancement. In self-attention, G acts as both a key and a value. Used as the query value.
7. The joint alignment and quantification method for large-scale mass spectrometry cohorts according to claim 1, characterized in that, The joint loss function simultaneously optimizes all or part of the parameters of the one-dimensional monotonic radial basis function space transformation network, the explicit intensity drift correction layer, the two-stream attention mechanism, and the quantitative regression module in the same computational graph through a single backpropagation. The joint loss function for: in, , , , , and These are, respectively, the quantitative regression error term, the retention time alignment consistency constraint, the quality control stability and detrending constraint, the time transformation monotonicity regularization term, the time transformation smoothness regularization term, and the L2 norm regularization constraint on the network parameters. , , , , and All are non-negative weight coefficients; the weight coefficients of each item in the joint loss function The method employs uncertainty-weighted, adaptive balancing, or gradient norm-based equalization to dynamically update the data during training. Time transformation function Alignment time grid after discretization The first derivative, Time transformation function Alignment time grid after discretization The second derivative; The quantitative regression error term includes label-based supervised loss and self-supervised consistency loss; The quality control stability and detrending constraints include one or a combination of the following: statistical dispersion minimization constraint, temporal correlation decoupling constraint, and adaptive baseline constraint; The statistical dispersion minimization constraint is as follows: for a preset target substance or feature set, calculate at least one of the following dispersion measures among the variance, standard deviation, coefficient of variation, and mean absolute deviation of the intensity of the quality control sample after correction, and minimize the dispersion measure to achieve response stability across time and / or across batches. The time-series correlation decoupling constraint is to calculate at least one of the following correlation measures: Pearson correlation coefficient, Spearman correlation coefficient, or mutual information, which measures the intensity of quality control samples as a function of injection order, collection time, or operation number, and to minimize the absolute value, square, or weighted term of the correlation measure in order to reduce the time-related system drift trend. The adaptive baseline constraint is as follows: when there is a lack of external quality control samples or the quality control samples are significantly different from the biological sample matrix, making the quality control constraint inapplicable, a robust statistical baseline based on the biological sample population is constructed, and the constraint objective is configured to minimize the systematic offset of the biological samples relative to the statistical baseline. The robust statistical baseline includes at least one of the sliding window median, quantile baseline, or low-rank approximation background.
8. The joint alignment and quantification method for large-scale mass spectrometry cohorts according to any one of claims 1-7, characterized in that, The joint alignment and quantification method employs a memory-memory hierarchical management and end-to-end sparse computation strategy, which includes: Streaming micro-batch loading: Employing a retention time sliding window mechanism and sample micro-batch strategy, only a portion of the sample data required for the current calculation that is within the retention time sliding window is loaded into memory, instead of loading the complete sample data all at once. The processing of the complete sample data is completed through iterative looping. Sparse compression and dynamic mapping: Micro-batch tensors are stored using coordinate lists or compressed sparse row formats, and high-dimensional sparse physical m / z coordinates are mapped to compact logical feature channel indices using dynamic ROI mapping technology, hash mapping or index lookup tables, thereby filtering out more than 99% of invalid background zero values before entering the computation graph. Sparse-aware computation: Sparse-aware computation is performed in the end-to-end deep learning network to construct a binary mask corresponding to non-zero features. The mask is used to shield invalid regions in the Softmax step of the Transformer self-attention mechanism, and the mask-filling mechanism is used to block gradient updates of zero-value regions during backpropagation, so as to focus on the chemical feature learning of effective ions. Adaptive resource scheduling: Real-time monitoring of GPU memory usage. When memory usage exceeds a preset threshold, a dynamic adjustment mechanism is automatically triggered to perform at least one of the following operations: reducing the micro-batch size, shortening the length of the retention time sliding window, or dynamically increasing the sparsity noise threshold.
9. A computer device, comprising: A memory, a processor, and a computer program stored in the memory and executable on the processor, characterized in that the processor executes the computer program to implement the joint alignment and quantification method for large-scale mass spectrometry queues as described in any one of claims 1-8.
10. A computer program product, comprising a computer program, characterized in that, When executed by a processor, the computer program implements the joint alignment and quantification method for large-scale mass spectrometry queues as described in any one of claims 1-8.