A method for batch effect subtraction of gene alternatively spliced isoforms

Abstract

The embodiment of the present disclosure provides a gene variable splicing isomer batch effect deduction method, which belongs to the technical field of data identification and specifically comprises the following steps: acquiring a GTEx transcript count matrix as a parameter estimation reference data source; estimating a negative binomial distribution parameter and sampling to generate a new count matrix; injecting a batch effect and a differential transcript to use a DTU signal to construct a simulation data set; performing an ILR space countermeasure batch process on the simulation data set, that is, constructing an intragenic composition vector with a gene as a unit and performing ILR transformation, learning a batch discrimination direction in the ILR space with a batch variable as a supervision signal and projecting to remove, and then recovering to a count matrix after batch removal through inverse ILR transformation; and running a differential analysis method before and after batch removal to obtain a result table and comparing and evaluating. Through the scheme of the present disclosure, false positives are effectively reduced, the accuracy and robustness of DTU detection are improved, and the scheme is suitable for transcriptome data analysis.

Description

Technical Field

[0001] This disclosure relates to the field of bioinformatics technology, and in particular to a method for batch effect subtraction of gene alternative splicing isoforms. Background Technology

[0002] Currently, batch effect is one of the main factors affecting the reliability of statistical tests and the reproducibility of results in differential analysis of transcriptome data. Batch effect arises from non-biological differences such as library construction, sequencing platform and process, reagent batch number, sequencing depth, and sample processing conditions, which can introduce systematic bias and lead to false positives or false negatives in differential analysis results.

[0003] The most common de-batch strategy currently used is to introduce batch covariates into the difference analysis model by designing a matrix, and then subtract the batch effect during the fitting and testing process. This approach is widely used in various difference analysis frameworks, such as satuRn, limma diffsplice, and DEXSeq.

[0004] However, in differential analysis at the isoform level, many methods use the relative usage ratio of transcripts within the same gene as the core analytical object. This type of data has compositional constraints (such as a ratio sum of 1) and exhibits strong coupling correlations. When batches are directly subtracted through a design matrix, the correction process may alter the overall proportional structure within the gene, thereby mistaking the technical effect of the perturbed proportional structure as a true differential signal, thus increasing the risk of false positives, especially when batch effects are strong or correlated with grouping.

[0005] It is evident that there is an urgent need for a batch effect subtraction method for gene alternative splicing isoforms that can reduce the impact of batch effects while minimizing the disruption of the proportional structure within the gene. Summary of the Invention

[0006] In view of this, the present disclosure provides a method for batch effect subtraction of gene alternative splicing isoforms, which at least partially solves some of the problems existing in the prior art.

[0007] This disclosure provides a method for batch effect subtraction of gene alternative splicing isoforms, including: Step 1: Obtain the transcript count matrix, sample information table, and tx2gene mapping table for the entire tissue and sample of GTEx, and use the transcript count matrix as the reference data source for parameter estimation. Step 2: Calculate the count distribution parameters of the reference count matrix by transcript dimension, including the mean parameter of each transcript, and further estimate the discrete parameters of the negative binomial distribution; Step 3: Perform negative binomial random sampling on each transcript based on the negative binomial distribution parameter set to generate a new transcript counting matrix; Step 4: Inject batch effect and DTU signals into the new transcript count matrix according to preset rules to obtain a simulated dataset containing batch factors and true signals. Step 5: Align the count matrix and sample information table in the simulated dataset with sample keys, align the count matrix and tx2gene with transcript keys, and complete the consistency check to obtain the aligned input dataset; Step 6: Determine the differential analysis grouping variable "condition" and the batch variable "batch" from the sample information table in the input dataset to form a modeling information table containing sample key, condition, and batch. Step 7: Filter the centrally aligned count matrix of the input dataset to obtain the original input matrix for difference analysis; Step 8: Based on the modeling information table and the original input matrix, perform the ILR spatial adversarial debaterization process to generate the debaterized input matrix; Step 9: Run the preset difference analysis method on the original input matrix to obtain the baseline results table; Step 10: Run the target difference analysis method on the batch input matrix to obtain the analysis results table; Step 11: Summarize and output the baseline results table and analysis results table. The output should include at least the results tables for each method and evaluation / visualization files for comparison.

[0008] According to a specific implementation of an embodiment of this disclosure, step 1 specifically includes: The transcript count matrix of the entire tissue and sample from GTEx was used as the reference data source for negative binomial parameter estimation. A certain number of sample columns were extracted from it for subsequent estimation of the mean and discrete parameters of each transcript. At the same time, the mapping table tx2gene between transcripts and genes was obtained.

[0009] According to a specific implementation of an embodiment of this disclosure, step 4 specifically includes: Step 4.1: Apply a systematic shift to the transcript composition structure within a gene, taking genes as the unit. Between preset batch groups, proportionally shift the relative usage ratio of transcripts of some genes to obtain a count matrix with batch confounding and simultaneously generate corresponding batch labels. Step 4.2: On the batch counting matrix, select the genes to be injected with DTU according to a preset ratio, and replace or redistribute the usage structure of several transcripts within each selected gene to form a clear differential transcript usage signal between different conditions, while recording the injected genes and transcript sets.

[0010] According to a specific implementation of an embodiment of this disclosure, step 5 specifically includes: Step 5.1: Align the counting matrix and the sample information table according to the sample key to ensure that each sample column of the counting matrix has a corresponding record in the sample information table; Step 5.2: Align the transcript rows of the counting matrix with the tx2gene mapping table by transcript key to ensure that each transcript can be mapped to the corresponding gene, and process missing, duplicate or unmatched entries to obtain the aligned input dataset.

[0011] According to a specific implementation of an embodiment of this disclosure, step 7 specifically includes: Low-expression transcript filtering and intragene transcript count constraints were performed on the aligned count matrix to obtain the original input matrix for differential analysis, while keeping its sample columns consistent with the modeling information table.

[0012] According to a specific implementation of an embodiment of this disclosure, step 8 specifically includes: Step 8.1: Construct an intragene composition vector by counting transcripts within the same gene, using genes as the unit, to obtain an intragene composition representation; Step 8.2: Perform ILR transformation on the gene composition vector to obtain the ILR space representation matrix; Step 8.3: Using the batch variable in the modeling information table as the supervision signal, the batch discrimination direction of matrix learning / estimation is represented in the ILR space; Step 8.4: Project and remove the components of the ILR space representation matrix in the batch discrimination direction to obtain the de-batch ILR space representation matrix; Step 8.5: Perform an inverse ILR transformation on the debatable ILR spatial representation to obtain the debatable compositional representation and restore it to the input form usable for difference analysis, forming the debatable input matrix.

[0013] According to a specific implementation of an embodiment of this disclosure, step 8.2 specifically includes: For any gene's intragene composition vector Taking the element-wise logarithm yields Then its ILR space representation matrix Written as:

[0014] Where H is a predetermined orthogonal basis transformation matrix.

[0015] According to a specific implementation of an embodiment of this disclosure, step 8.3 specifically includes: The representation of each sample in the ILR space is used as a feature vector, and the batch value corresponding to that sample is used as a supervision label to form a dataset for batch discrimination training. The discrimination direction is obtained by maximizing inter-batch differences and minimizing intra-batch differences, wherein the expression for the discrimination direction is:

[0016] in, Represents the batch-to-batch scatter matrix. Represents the scatter matrix within the batch. To determine the direction of the learned batch, This indicates the transpose operation.

[0017] According to a specific implementation of an embodiment of this disclosure, step 8.4 specifically includes: Based on the learned batch discrimination direction, the components of the ILR space representation matrix in that direction are projected and removed to obtain the de-batched ILR space representation matrix: .

[0018] According to a specific implementation of an embodiment of this disclosure, step 8.5 specifically includes: Let z be the debatable ILR spatial representation vector. First, perform an inverse ILR transformation to restore it to the intragene composition space, then scale it according to the total gene count to restore it to the transcript count form, thus forming the debatable input matrix. The expression for the inverse ILR transformation is:

[0019] in, This indicates element-wise exponentiation. The matrix used for the inverse ILR transformation. Let be the normalization coefficient, so that The sum of all components is 1.

[0020] The batch effect subtraction scheme for gene alternative splicing isoforms in this embodiment includes: Step 1, obtaining the transcript count matrix, sample information table, and tx2gene mapping table of the entire tissue and sample from GTEx, and using the transcript count matrix as the reference data source for parameter estimation; Step 2, calculating the count distribution parameters of the reference count matrix according to the transcript dimension, including the mean parameter of each transcript, and further estimating the discrete parameters of the negative binomial distribution; Step 3, performing negative binomial random sampling on each transcript based on the negative binomial distribution parameter set to generate a new transcript count matrix; Step 4, injecting batch effect and DTU signals into the new transcript count matrix according to preset rules to obtain a simulated dataset containing batch factors and true signals; Step 5, aligning the count matrix and sample information table in the simulated dataset with sample keys, aligning the count matrix and tx2gene with transcript keys, and completing a... Step 6: Perform a consistency check to obtain the aligned input dataset; Step 7: Determine the differential analysis grouping variable "condition" and batch variable "batch" from the sample information table in the input dataset, forming a modeling information table containing sample keys, condition, and batch; Step 8: Filter the aligned count matrix in the input dataset to obtain the original input matrix for differential analysis; Step 9: Based on the modeling information table and the original input matrix, perform an ILR spatial adversarial debaterization process to generate a debaterized input matrix; Step 10: Run a preset differential analysis method on the original input matrix to obtain a baseline result table; Step 11: Run a target differential analysis method on the debaterized input matrix to obtain an analysis result table; Step 22: Summarize and output the baseline result table and analysis result table, with the output including at least the result tables for each method and evaluation / visualization files for comparison.

[0021] The beneficial effects of the embodiments of this disclosure are as follows: Through the innovative combination of "ILR transformation + adversarial batch de-batching" in the scheme of this disclosure, the batch effect is efficiently deducted while protecting the compositional structure of the gene, which effectively reduces false positives, improves the detection efficiency of differential transcript usage analysis, and shows good robustness and scalability under different sample sizes. Attached Figure Description

[0022] To more clearly illustrate the technical solutions of the embodiments of this disclosure, the accompanying drawings used in the embodiments will be briefly introduced below. Obviously, the accompanying drawings described below are only some embodiments of this disclosure. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0023] Figure 1 A flowchart illustrating a batch effect subtraction method for gene alternative splicing isoforms provided in this embodiment of the disclosure; Figure 2 This is a schematic diagram illustrating the specific implementation process of a batch effect subtraction method for gene alternative splicing isoforms provided in this embodiment of the disclosure. Detailed Implementation

[0024] The embodiments of this disclosure will now be described in detail with reference to the accompanying drawings.

[0025] The following specific examples illustrate the implementation of this disclosure. Those skilled in the art can easily understand other advantages and effects of this disclosure from the content disclosed in this specification. Obviously, the described embodiments are only a part of the embodiments of this disclosure, and not all of them. This disclosure can also be implemented or applied through other different specific embodiments, and the details in this specification can also be modified or changed based on different viewpoints and applications without departing from the spirit of this disclosure. It should be noted that, in the absence of conflict, the following embodiments and features in the embodiments can be combined with each other. Based on the embodiments in this disclosure, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of this disclosure.

[0026] It should be noted that various aspects of embodiments within the scope of the appended claims are described below. It will be apparent that the aspects described herein can be embodied in a wide variety of forms, and any particular structure and / or function described herein is merely illustrative. Based on this disclosure, those skilled in the art will understand that one aspect described herein can be implemented independently of any other aspect, and two or more of these aspects can be combined in various ways. For example, any number of aspects set forth herein can be used to implement the device and / or practice the method. Additionally, this device and / or method can be implemented using structures and / or functionalities other than one or more of the aspects set forth herein.

[0027] It should also be noted that the illustrations provided in the following embodiments are only schematic representations of the basic concept of this disclosure. The illustrations only show the components related to this disclosure and are not drawn according to the number, shape and size of the components in actual implementation. In actual implementation, the form, quantity and proportion of each component can be arbitrarily changed, and the layout of the components may also be more complex.

[0028] Furthermore, specific details are provided in the following description to facilitate a thorough understanding of the examples. However, those skilled in the art will understand that the described aspects can be practiced without these specific details.

[0029] This disclosure provides a method for batch effect subtraction of gene alternative splicing isoforms, which can be applied to transcriptome data differential analysis in bioinformatics engineering scenarios.

[0030] See Figure 1 This is a flowchart illustrating a batch effect subtraction method for gene alternative splicing isoforms provided in an embodiment of this disclosure. Figure 1 As shown, the method mainly includes the following steps: Step 1: Obtain the transcript count matrix, sample information table, and tx2gene mapping table for the entire tissue and sample of GTEx, and use the transcript count matrix as the reference data source for parameter estimation. In practice, the data acquisition process is as follows: the transcript count matrix of the entire tissue and sample of GTEx is used as the reference data source (parameter pool) for negative binomial (NB) parameter estimation. A certain number of sample columns are extracted from it for subsequent estimation of the mean and discrete parameters of each transcript. At the same time, the transcript-gene mapping table tx2gene is obtained (or loaded) to determine the transcript composition structure within the gene and to serve as the structural basis for the subsequent simulation count matrix construction and injection steps.

[0031] Step 2: Calculate the count distribution parameters of the reference count matrix by transcript dimension, including the mean parameter of each transcript, and further estimate the discrete parameters of the negative binomial distribution; In practice, based on the GTEx reference counting matrix obtained in step 1, the mean parameter is calculated according to the transcript dimension, and the discrete parameters of the negative binomial distribution are estimated by combining statistics such as variance, forming the NB parameter set (parameter pool) for subsequent NB random sampling to generate the simulated counting matrix.

[0032] Step 3: Perform negative binomial random sampling on each transcript based on the negative binomial distribution parameter set to generate a new transcript counting matrix; In practice, based on the set of negative binomial distribution parameters obtained in step 2, each transcript is randomly sampled according to its corresponding NB parameter to generate a new transcript counting matrix. This simulated counting matrix is ​​consistent with the GTEx reference data in terms of the overall counting distribution, but the sample level is obtained by resampling, thus serving as the basic input for subsequent batch effects and DTU signal injection.

[0033] Step 4: Inject batch effect and DTU signals into the new transcript count matrix according to preset rules to obtain a simulated dataset containing batch factors and true signals. In practice, the batch effect and DTU signal injection are performed in the following manner: A. Batch effect injection: Apply a systematic shift to the transcript composition structure within a gene, on a gene-by-gene basis. That is, between preset batches, the transcripts of some genes are proportionally shifted according to their relative usage ratio. The shift intensity is controlled by preset parameters (e.g., whether the batch direction is related to the condition, the shift amplitude, and the ratio of genes / transcripts involved in the shift), thereby obtaining a count matrix with batch contamination and simultaneously generating corresponding batch labels.

[0034] B. DTU signal injection: On the batch counting matrix, select genes to be injected with DTU according to a preset ratio, and replace / redistribute the usage structure of several transcripts within each selected gene (e.g., exchange the relative proportions of some transcripts within the gene) to form a clear differential transcript usage signal between different conditions. At the same time, record the injected gene and transcript set as ground truth labeling information for subsequent evaluation.

[0035] Step 5: Align the count matrix and sample information table in the simulated dataset with sample keys, align the count matrix and tx2gene with transcript keys, and complete the consistency check to obtain the aligned input dataset; In practice, the counting matrix obtained in step 4 is aligned with the sample information table by the sample key to ensure that each sample column of the counting matrix has a corresponding record in the sample information table; at the same time, the transcript row of the counting matrix is ​​aligned with the tx2gene mapping table by the transcript key to ensure that each transcript can be mapped to the corresponding gene; and missing, duplicate or unmatched entries are processed to obtain an aligned input dataset that can be used for subsequent analysis.

[0036] Step 6: Determine the differential analysis grouping variable "condition" and the batch variable "batch" from the sample information table in the input dataset to form a modeling information table containing sample key, condition, and batch. In practice, based on the aligned sample information table in the input dataset obtained in step 5, the difference analysis grouping variable condition and batch variable batch are determined, and a modeling information table containing sample key, condition, and batch is extracted and formed for subsequent design matrix construction and method invocation.

[0037] Step 7: Filter the centrally aligned count matrix of the input dataset to obtain the original input matrix for difference analysis; In practice, necessary expression filtering and structure filtering (such as low-expression transcript filtering and intragenic transcript quantity constraints) are performed on the centrally aligned count matrix of the input dataset obtained in step 5 to obtain the original input matrix (without batch removal) for differential analysis, and its sample columns are kept consistent with the modeling information table in step 6.

[0038] Step 8: Based on the modeling information table and the original input matrix, perform the ILR spatial adversarial debaterization process to generate the debaterized input matrix; In practice, the ILR spatial adversarial debaterization process is executed, and the debaterization input matrix is ​​generated. This can be achieved through the following steps: (1) Constructing and normalizing the gene composition: Using genes as the unit, the transcript count of the same gene is converted into the gene composition ratio in each sample, which is used for subsequent stable transformation and batch processing in the composition space. The specific process is as follows: Based on the original input matrix obtained in step 7, the transcript counts are grouped by gene according to the tx2gene mapping relationship. Within each sample, the transcript count vector of the same gene is normalized to the gene composition vector (i.e., each transcript count is divided by the total count of the gene in the sample), thus obtaining the gene composition representation for subsequent ILR transformation and de-batch processing. (2) ILR transformation to obtain ILR space representation: Perform ILR transformation on the gene composition vector obtained in step (1) to map it from the composition space with sum and constraint to the Euclidean space without sum and constraint, and obtain the ILR space representation matrix of each sample on each gene. The specific process is as follows: Based on the obtained intragenic composition representation, an ILR transformation is performed on the composition vector of each gene, mapping it from a composition space with sum and constraint to an unconstrained Euclidean space, resulting in the corresponding ILR space representation matrix, which is used for subsequent batch discrimination direction learning and projection removal. For any intragenic composition vector denoted as x (where each component is positive and sums to 1), the element-wise logarithm is taken to obtain... Then its ILR space representation It can be written as:

[0039] Where H is a predetermined orthogonal basis transformation matrix; (3) Adversarial learning of batch discrimination directions: Using the batch variable in the sample information table as the supervision signal, a batch discrimination model is fitted on the ILR space representation matrix to obtain the discrimination directions that can distinguish different batches to the greatest extent (when the batch is multi-level, the corresponding set of discrimination directions or its spanned subspace is obtained). The discrimination directions are used to characterize the main batch components in the ILR space. The specific process is as follows: Using the batch variable determined in step 6 as the supervision signal, and the obtained ILR space representation as the input feature, a batch discrimination model is fitted to learn the discrimination direction that can distinguish different batches to the greatest extent. The discrimination direction is used to characterize the main batch components in the ILR space and serves as the basis for subsequent projection removal. Specifically, it can be implemented in the following way: (1) Sample and feature organization: The representation of each sample in the ILR space is used as a feature vector, and the batch value corresponding to the sample is used as a supervision label to form a dataset for batch discrimination training.

[0040] (2) Determining the target and direction: The target and direction are determined by maximizing the inter-batch difference and minimizing the intra-batch difference. The target and direction can be expressed as:

[0041] in, Represents the batch-to-batch scatter matrix. Represents the scatter matrix within the batch. To determine the direction of the learned batch, Indicates the transpose operation; (4) Projection to remove batch components: Project and remove the components of the ILR space representation obtained in step (2) on the batch discrimination direction (or its spanned subspace) obtained in step (3) to obtain the batch-removed ILR space representation matrix. The specific process is as follows: Based on the learned batch discrimination direction, the component of the obtained ILR space representation in that direction is projected and removed, thereby weakening the main batch-driven change component in the ILR space, resulting in the de-batch ILR space representation matrix. For any sample's ILR feature vector z, its de-batch representation can be obtained as follows: first, calculate its projection component in direction w; then, the representation after removing this direction component is... ; For cases with a single discrimination direction, the component of the ILR vector of each sample in that direction can be removed to achieve orthogonal projection removal; (5) Inverse ILR transformation and restoration to a usable input matrix: Let z be the debatable ILR spatial representation from step (4). Perform inverse ILR transformation to restore the debatable intragene composition representation. Then, using the total gene count of each sample as the scale, restore the composition representation to the transcript count form to obtain the debatable input matrix that can be directly used for subsequent differential analysis. The expression for the inverse ILR transformation is:

[0042] in, This indicates element-wise exponentiation. The matrix used for the inverse ILR transformation. Let be the normalization coefficient, so that The sum of all components is 1.

[0043] Step 9: Run the preset difference analysis method on the original input matrix to obtain the baseline results table; In practice, differential analysis can be run on the original input matrix (without batch removal) obtained in step 7 to obtain a baseline result table; the baseline result table contains at least the p-value for each transcript / event, and further performs multiple test corrections to obtain the FDR.

[0044] Step 10: Run the target difference analysis method on the batch input matrix to obtain the analysis results table; In practice, the same differential analysis method can be run on the batch-free input matrix obtained in step 8 to obtain the result table of the present invention; the result table of the present invention includes at least the p-value of each transcript / event, and further performs multiple test correction to obtain FDR.

[0045] Step 11: Summarize and output the baseline results table and analysis results table. The output should include at least the results tables for each method and evaluation / visualization files for comparison.

[0046] In practice, the two sets of results obtained in steps 9 and 10 are summarized and output in a unified manner. The output content includes at least the baseline result table and the result table of this invention (including fields such as p-value and FDR), as well as evaluation and visualization files for comparative analysis.

[0047] The batch effect subtraction method for gene alternative splicing isoforms provided in this embodiment introduces ILR transformation to process compositional data: Addressing the constraint that the proportion of transcripts within a gene has a "sum of 1" characteristic, an isometric logarithmic ratio (ILR) transformation is used to map the data from the compositional space to an unconstrained Euclidean space. This transformation preserves the relative proportions between transcripts, laying the foundation for accurate batch effect subtraction in Euclidean space, thus avoiding the false positive problem caused by the destruction of the intrageneal proportional structure due to direct batch subtraction in traditional methods. An adversarial learning strategy is employed in the ILR space to subtract batch effects: using batch variables as supervisory signals, the method learns the discrimination direction that best distinguishes different batches in the ILR space, and then removes the component of the data in that direction through orthogonal projection. This adversarial mechanism can accurately identify and isolate batch-related variations, effectively preserving the true differential transcript usage signal even when batch effects are strong or confounded with biological groupings, significantly improving the sensitivity and accuracy of detection. A complete simulation validation framework is constructed: parameters are estimated based on GTEx real data, and batch effects and differential signals are controllably injected to generate a benchmark dataset with ground truth labels. This framework provides objective standards for evaluating the performance of methods, ensuring their reliability and reproducibility in practical applications.

[0048] The method of this disclosure will be further described below with reference to a specific embodiment. The invention will be evaluated on semi-analog data GTX. Tools that are compared with this invention include satuRn, DRIMSeq, DEXseq, limma_diffsplice, and edgR_diffsplice. 1) Experimental analysis and prediction were conducted on 40 samples, 20 case samples, and 20 control samples; 2) Experimental analysis and prediction were performed on 100 samples, 50 case samples and 50 control samples.

[0049] In the two experiments described above, under the same edgeR filtering and count input conditions, when comparing benchmark datasets with sample sizes of 20 cases / 20 control and 50 cases / 50 control, AC_satuRn consistently showed higher TPR across all FDR threshold ranges, especially in the low FDR region (near commonly used thresholds such as 0.01, 0.05, and 0.1), demonstrating a more significant advantage. This indicates that the method of this invention can identify more true DTU signals while strictly controlling for false positives. As the sample size increased from 20 vs 20 to 50 vs 50, the overall TPR of each method increased, but the AC_satuRn curve remained consistently above the main control method, indicating that the method has good robustness and scalability across different sample sizes. In contrast, some control methods showed limited TPR improvement in the low FDR region, reflecting limited detection sensitivity in the presence of batch contamination.

[0050] As shown in Table 1, under both 20 case / 20 control and 50 case / 50 control sample size settings, AC_satuRn achieved the highest value in the area under the FDR-TPR curve (AUC_mean), indicating that the method of this invention has superior overall detection performance when considering different FDR thresholds. With the increase in sample size, the AUC_mean of all methods improved, but the improvement in AC_satuRn was more significant, further demonstrating that it can more stably retain the true DTU signal after removing batch effects, while avoiding performance loss caused by batch interference in proportional modeling.

[0051] Table 1

[0052] It should be understood that the various parts of this disclosure can be implemented in hardware, software, firmware, or a combination thereof.

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

Claims

1. A method for batch effect subtraction of gene alternative splicing isoforms, characterized in that, include: Step 1: Obtain the transcript count matrix, sample information table, and tx2gene mapping table for the entire tissue and sample of GTEx, and use the transcript count matrix as the reference data source for parameter estimation. Step 2: Calculate the count distribution parameters of the reference count matrix by transcript dimension, including the mean parameter of each transcript, and further estimate the discrete parameters of the negative binomial distribution; Step 3: Perform negative binomial random sampling on each transcript based on the negative binomial distribution parameter set to generate a new transcript counting matrix; Step 4: Inject batch effect and DTU signals into the new transcript count matrix according to preset rules to obtain a simulated dataset containing batch factors and true signals. Step 5: Align the count matrix and sample information table in the simulated dataset with sample keys, align the count matrix and tx2gene with transcript keys, and complete the consistency check to obtain the aligned input dataset; Step 6: Determine the differential analysis grouping variable "condition" and the batch variable "batch" from the sample information table in the input dataset to form a modeling information table containing sample key, condition, and batch. Step 7: Filter the centrally aligned count matrix of the input dataset to obtain the original input matrix for difference analysis; Step 8: Based on the modeling information table and the original input matrix, perform the ILR spatial adversarial debaterization process to generate the debaterized input matrix; Step 9: Run the preset difference analysis method on the original input matrix to obtain the baseline results table; Step 10: Run the target difference analysis method on the batch input matrix to obtain the analysis results table; Step 11: Summarize and output the baseline results table and analysis results table. The output should include at least the results tables for each method and evaluation / visualization files for comparison.

2. The method according to claim 1, characterized in that, Step 1 specifically includes: The transcript count matrix of the entire tissue and sample from GTEx was used as the reference data source for negative binomial parameter estimation. A certain number of sample columns were extracted from it for subsequent estimation of the mean and discrete parameters of each transcript. At the same time, the mapping table tx2gene between transcripts and genes was obtained.

3. The method according to claim 1, characterized in that, Step 4 specifically includes: Step 4.1: Apply a systematic shift to the transcript composition structure within a gene, taking genes as the unit. Between preset batch groups, proportionally shift the relative usage ratio of transcripts of some genes to obtain a count matrix with batch confounding and simultaneously generate corresponding batch labels. Step 4.2: On the batch counting matrix, select the genes to be injected with DTU according to a preset ratio, and replace or redistribute the usage structure of several transcripts within each selected gene to form a clear differential transcript usage signal between different conditions, while recording the injected genes and transcript sets.

4. The method according to claim 1, characterized in that, Step 5 specifically includes: Step 5.1: Align the counting matrix and the sample information table according to the sample key to ensure that each sample column of the counting matrix has a corresponding record in the sample information table; Step 5.2: Align the transcript rows of the counting matrix with the tx2gene mapping table by transcript key to ensure that each transcript can be mapped to the corresponding gene, and process missing, duplicate or unmatched entries to obtain the aligned input dataset.

5. The method according to claim 1, characterized in that, Step 7 specifically includes: Low-expression transcript filtering and intragene transcript count constraints were performed on the aligned count matrix to obtain the original input matrix for differential analysis, while keeping its sample columns consistent with the modeling information table.

6. The method according to claim 1, characterized in that, Step 8 specifically includes: Step 8.1: Construct an intragene composition vector by counting transcripts within the same gene, using genes as the unit, to obtain an intragene composition representation; Step 8.2: Perform ILR transformation on the gene composition vector to obtain the ILR space representation matrix; Step 8.3: Using the batch variable in the modeling information table as the supervision signal, the batch discrimination direction of matrix learning / estimation is represented in the ILR space; Step 8.4: Project and remove the components of the ILR space representation matrix in the batch discrimination direction to obtain the de-batch ILR space representation matrix; Step 8.5: Perform an inverse ILR transformation on the debatable ILR spatial representation to obtain the debatable compositional representation and restore it to the input form usable for difference analysis, forming the debatable input matrix.

7. The method according to claim 6, characterized in that, Step 8.2 specifically includes: For any gene's intragene composition vector Taking the element-wise logarithm yields Then its ILR space representation matrix Written as: Where H is a predetermined orthogonal basis transformation matrix.

8. The method according to claim 6, characterized in that, Step 8.3 specifically includes: The representation of each sample in the ILR space is used as a feature vector, and the batch value corresponding to that sample is used as a supervision label to form a dataset for batch discrimination training. The discrimination direction is obtained by maximizing inter-batch differences and minimizing intra-batch differences, wherein the expression for the discrimination direction is: in, Represents the batch-to-batch scatter matrix. Represents the scatter matrix within the batch. To determine the direction of the learned batch, This indicates the transpose operation.

9. The method according to claim 6, characterized in that, Step 8.4 specifically includes: Based on the learned batch discrimination direction, the components of the ILR space representation matrix in that direction are projected and removed to obtain the de-batched ILR space representation matrix: 。 10. The method according to claim 6, characterized in that, Step 8.5 specifically includes: Let z be the debatable ILR spatial representation vector. First, perform an inverse ILR transformation to restore it to the intragene composition space, then scale it according to the total gene count to restore it to the transcript count form, thus forming the debatable input matrix. The expression for the inverse ILR transformation is: in, This indicates element-wise exponentiation. The matrix used for the inverse ILR transformation. Let be the normalization coefficient, so that The sum of all components is 1.

Images