A gene / locus screening method based on data discretization and heuristic search loop

Abstract

This invention discloses a gene / locus screening method based on data discretization and a heuristic search loop, aiming to screen for the fewest genes / locus combinations to achieve the greatest differentiation of sample types. The method compares the discretized data of the current loop with the previous round's best gene / locus combination, selects the optimal combination, and continuously adjusts hyperparameters based on the score results, iteratively searching for the gene / locus combination with the fewest genes / locus combinations and the highest optimal score. The algorithmic logic of this invention forces the optimization of the gene / locus screening process towards minimizing the number of genes / locus combinations and maximizing the differentiation of sample types, thereby making the screened gene / locus combinations more efficient and specific.

Description

Technical Field

[0001] This invention belongs to the field of bioinformatics technology, specifically relating to a computer algorithm based on data discretization and heuristic search loops that can be applied to biological gene / site screening logic. It is suitable for screening gene / site combinations and ultimately applied to disease screening and classification. Background Technology

[0002] Gene screening is the process of selecting genes from a large pool of genes that are associated with specific biological processes, diseases, or phenotypes through experimental or computational analysis. The goal of gene screening is to identify genes that play a crucial role in disease development, cellular function regulation, or other biological phenomena, providing candidate genes for subsequent functional validation or clinical applications.

[0003] Gene screening is mainly divided into two categories:

[0004] Experimental screening: This involves directly intervening in gene function through experimental methods, such as RNA interference (RNAi) or CRISPR-Cas9 gene editing technology, to induce loss of function or enhancement in thousands of genes. By observing phenotypic changes in cells or organisms, researchers can determine which genes play important roles under specific conditions. This type of method is often used to discover potential therapeutic targets related to cancer, neurodegenerative diseases, and other conditions.

[0005] Computational screening: This involves analyzing gene expression data, biomarkers, or other high-throughput omics data using statistical and machine learning methods. Techniques such as differential gene expression analysis, association analysis, and cluster analysis can be used to screen for genes significantly associated with diseases or specific phenotypes from thousands of genes. This type of screening is often used to extract information from large-scale datasets, such as screening for key driver genes in cancer using RNA sequencing data.

[0006] In gene screening, the importance of combinatorial gene biomarkers is becoming increasingly prominent. These biomarkers are not simply a collection of genes associated with a specific disease, but rather, through precise screening and optimization, combinations of genes that can reflect the occurrence, development, and individual differences of diseases with higher specificity and accuracy. For example, compared to single gene biomarkers, combinatorial gene biomarkers can significantly reduce false positive and false negative rates and improve the diagnostic and predictive capabilities for complex diseases. This approach is particularly suitable for multifactorial driven diseases, such as cancer and autoimmune diseases, enabling a more comprehensive capture of their underlying molecular mechanisms, thus laying a solid foundation for personalized treatment and precision medicine.

[0007] Commonly used methods for computational gene screening include differential expression analysis, cluster analysis, machine learning, and gene network analysis. For example, differential expression analysis uses statistical tools (such as DESeq2) to screen for genes that show significant changes in disease samples, but it is easily limited by sample size and false positive rate; cluster analysis (such as K-means) is used to discover co-regulatory patterns of genes, but it is sensitive to parameters and susceptible to noise; machine learning methods such as random forests can screen key genes based on feature importance, but they are prone to overfitting in small sample data; gene network analysis screens core genes by constructing co-expression networks, but it relies on prior knowledge of network topology. Although these gene screening methods can screen for disease-related genes, they generally suffer from drawbacks such as high false positive rate, low specificity, sensitivity to parameters or sample size, limited number of sample types that can be distinguished, susceptibility to overfitting, and poor ability to capture nonlinear relationships.

[0008] In the biomedical field, there is often a discrepancy between theoretical predictions and actual experimental results. The same gene may be associated with multiple diseases, and individual patient differences exist. Therefore, designing reasonable computational methods to scientifically address these issues remains a challenge in current research. Summary of the Invention

[0009] This invention addresses the shortcomings of existing technologies by providing a computational gene / locus combination screening method. It efficiently achieves gene / locus screening through a gene / locus screening algorithm based on a heuristic search cycle, incorporating biologically required conditions into the screening process to achieve efficient, highly specific, and controllable gene / locus combination screening.

[0010] This invention discloses a gene / locus screening method based on data discretization and a heuristic search loop, aiming to screen for the fewest genes / locus combinations to achieve the greatest differentiation among sample types. The method uses data discretization and a scoring system as its core for gene / locus screening. Specifically, it compares the scoring results of the discretized data of the current loop with the previous round's best gene / locus combination, selecting the optimal combination and iteratively searching for the best combination. Data discretization involves normalizing the original data, specifically dividing the normalized data into multiple intervals, and discretizing the data belonging to each interval into a unique label. The scoring system concatenates the discretized data of a gene / locus combination according to sample type, obtaining concatenated discretized labels for all genes / locus combinations of each sample type. The number of unique concatenated discretized labels under the current gene / locus combination is analyzed, which is the scoring criterion of the scoring system. The core of the heuristic lies in the process of continuously adjusting hyperparameters and the final result based on the results of the scoring system. The hyperparameters include the discretization threshold, the number of genes / sites in the current combination, the gene / site number threshold, the number of no-improvement cycles, the no-improvement cycle threshold, and the number of updated / replaced genes / sites. The final result of the method is the gene / site combination with the fewest gene / sites and the largest optimal score.

[0011] The specific technical solution of this invention is as follows:

[0012] A gene / site screening method based on data discretization and heuristic search loop includes the following steps:

[0013] (1) Initialize and define hyperparameters

[0014] Discretization threshold: Discretization threshold x1~x w 0 <x1<……<x w <1, Discretization threshold x1~x w Divide the closed real number interval from 0 to 1 into discontinuous regions [0, x1], ..., [x...]. w [1], w represents even numbers from 2 to 10;

[0015] k: The number of genes / locus sites in the current combination; the value of k can be set manually as needed. The minimum value of k is the floor of log(Dis_C)Sam_C, and the maximum value is km. It is usually initialized to the floor of log(Dis_C)Sam_C, where Dis_C is the number of discretization intervals and Sam_C is the number of sample types to be distinguished.

[0016] km: The threshold for the number of genes / locus sites, which is the maximum number of genes / locus sites that are to be distinguished from the sample types. The value of km can be set manually as needed, usually 2 to 4 times the rounded-up value of log(Dis_C)Sam_C, and its minimum value is the rounded-up value of log(Dis_C)Sam_C, where Dis_C is the number of discretization intervals and Sam_C is the number of sample types to be distinguished.

[0017] m: No improvement cycle count, which is the number of consecutive cycles when there is no improvement during the scoring process. The initial value is defined as 1, and the maximum value is ct, which is the no improvement cycle threshold.

[0018] ct: No Improvement Cycle Threshold, which is the maximum number of consecutive cycles counted when there is no improvement during the scoring process; ct: The value can be set manually as needed, usually greater than or equal to X / 5, approximately 1000~80000, where X is the number of all selectable genes / sites.

[0019] n: The number of genes / sites to be replaced. The number of genes / sites in the optimal combination of genes / sites to be replaced. The initial value is defined as 1, and the maximum value is km-1.

[0020] (2) Discretization process:

[0021] Select several independent samples from each sample type in a bioinformatics database (which may be a publicly available bioinformatics database, a commercially available database that can be licensed, or a database constructed based on real independent samples). Obtain the raw data of all genes / locuses of all independent samples corresponding to each sample type. Normalize the raw data that is not in the range [0, 1] so that all data are scaled to within 0 to 1. Take the average of the data (raw data or normalized data within 0 to 1) of all independent samples corresponding to each gene / locus under each sample type (the raw data can be of the same type or different types). The sample type is selected from the subclass of the same category of features of disease, race, age, sex, or other biological or environmental factors that affect gene characteristics.

[0022] The different sample types in the method described in this invention are further subcategories under the same type of factors. For example, when the sample type is a disease, the different sample types described in this invention are different diseases (e.g., diabetes, hypertension, hyperlipidemia, etc., or lung cancer, liver cancer, stomach cancer, colorectal cancer, etc.); when the sample type is race, the different sample types described in this invention are different races (e.g., Caucasian, Asian, Brown, Black); when the sample type is age, the different sample types described in this invention are different age groups (e.g., children, youth, middle-aged, elderly); when the sample type is gender, the different sample types described in this invention are male or female; when the sample type is lifestyle or environmental factors, the different sample types described in this invention may include smokers and non-smokers, people with different exercise habits, or people with different dietary habits, etc.

[0023] The number of sample types can be set manually as needed, usually between 2 and 50.

[0024] The mean-normalized data is divided according to the initially defined discretization threshold. Each data point in the range [0, x1] is assigned a label a1, those in the range [x2, x3] are assigned a label a2, and so on, until a value in the range [x2, x3] is assigned a label a3. w ,1]] data is assigned the label a (w / 2+1) Data that does not fall within any of the above discretization thresholds are assigned empty labels. Each gene / locus has one label for each sample type. Genes / locus that have the same label for all sample types and all genes / locus containing empty labels are removed to obtain a standard gene / locus dataset.

[0025] (3) Scoring is based on a heuristic search loop algorithm, specifically including:

[0026] (i) Randomly select k genes / locus from the standard gene / locus dataset as the optimal gene / locus combination for initialization, and obtain the label spectrum formed by each label for each sample type corresponding to the optimal gene / locus combination. Compare the label spectrum of different sample types, calculate the number of sample types with non-repeating label spectrum, and define it as the optimal score.

[0027] (ii) Calculate the number of unique sample types in the label spectrum of the remaining k-1 gene / site combinations after removing any gene / site from the current best gene / site combination. This number is defined as the gene / site score of the removed gene / site. A total of k gene / site scores are obtained. The gene / site scores are sorted from largest to smallest. The top n genes with the highest rankings are selected to form the gene / site index to be replaced. From the standard dataset where the gene / sites in the best gene / site combination have been removed, n gene / sites are reselected to replace the gene / site index in the current best gene / site combination. Candidate gene / site combinations are constructed, and the label spectrum of each sample type corresponding to the candidate gene / site combination is obtained. The label spectra of different sample types are compared, and the number of unique sample types in the label spectrum is calculated. This number is defined as the candidate score.

[0028] (iii) Scoring Cycle:

[0029] If the candidate score is greater than the best score, the current candidate gene / site combination is replaced with the best gene / site combination, and the candidate score is replaced with the best score. Refer to step (ii) to construct a new candidate gene / site combination and calculate a new candidate score. Repeat the above operation until the updated candidate score is less than or equal to the best score of the previous round of replacement.

[0030] If the candidate score is less than or equal to the best score, keep the best gene / site combination unchanged, keep the best score unchanged, keep the gene / site score unchanged, keep the gene / site index to be replaced unchanged, and refer to step (ii) to construct a new candidate gene / site combination and calculate a new candidate score (i.e., reselect n genes / sites from the standard dataset that removes the genes / sites in the current best gene / site combination to replace the gene / site index, construct a candidate gene / site combination with the remaining kn gene combinations in the current best gene / site combination, and obtain the label spectrum of each sample type corresponding to the candidate gene / site combination formed by each label, compare the label spectrum of different sample types, calculate the number of sample types with non-repeating label spectrum, and define it as the candidate score), repeat the above operation, and update the m value. For each increase of the loop number, the m value is updated by 1 until the new candidate score is greater than the best score;

[0031] (4) When the scoring cycle in step (iii) reaches the point where the updated m value is greater than the ct value, the parameters are updated:

[0032] (a) If n < k - 1, save the current best gene / locus combination and the corresponding hyperparameters, reset the non - improvement cycle number m to 1, update the value of n by adding 1, keep the best gene / locus combination unchanged, keep the best score unchanged, keep the gene - locus scores unchanged, update the gene / locus index to be replaced according to the updated n value, construct a new candidate gene / locus combination and calculate the new candidate score according to step (ii), and enter the score loop (iii);

[0033] (b) If n = k - 1, adjust the discretization threshold, increase and / or decrease one or more values among x1 to x w and reset n to the initial value 1 and the non - improvement cycle number m to 1;

[0034] Repeat steps (2), (3) according to the new discretization threshold and update the parameters using steps (a) and (b) until the discretization threshold cannot be further adjusted, that is, when the latest discretization threshold divides the closed interval from 0 to 1 into completely continuous sub - intervals. If the current k is not greater than km, update the value of k by adding 1, reset the discretization threshold to the initial discretization threshold, reset n to 1, reset m to 1, and repeat steps (2), (3) and update the parameters using steps (a) and (b);

[0035] (5) Loop termination condition

[0036] When all loops and parameter updates reach the point where the current k is greater than km, exit all loops;

[0037] (6) Compare each best gene / locus combination saved in the loop and its corresponding hyperparameters, select the gene / locus combination with the largest best score value as the final target gene / locus combination. If there are multiple equal largest best score values, calculate the mean of the distances between adjacent regions of the closed real - number interval from 0 to 1 divided by the current discretization threshold for each candidate gene / locus combination, and select the gene / locus combination with the largest mean and / or the fewest genes / loci as the final target gene / locus combination.

[0038] In the method of the present invention, the original data is one or several of DNA sequencing data, RNA sequencing data, DNA methylation sequencing data, miRNA sequencing data, protein expression profile, DNA copy number variation profile, metabolomics data, genomics data, proteomics data, transcriptomics data, epigenomics data, immunomics data, pharmacomics data, systems biology data, phenomics data, lipidomics data, glycomics data, signalomics data, toxicomics data or evolutionary genomics data.

[0039] The raw data are selected from one or more of the following: methylation rate, gene expression level, mutation status, gene copy number, miRNA expression abundance, protein expression level, single nucleotide polymorphism frequency, structural variation copy number, epigenetic modification intensity, non-coding RNA expression abundance, phosphorylated protein ratio, protein-protein interaction intensity, metabolite concentration, metabolic pathway activity, immune cell infiltration ratio, cytokine concentration, microbial abundance, metabolite contribution, target inhibition efficiency, glucose metabolism concentration, or lipid metabolism concentration.

[0040] The method described in this invention is characterized in that the characteristics of the other biological or environmental factors affecting gene characteristics are selected from genetic background, tissue origin, lifestyle, and environmental exposure factors.

[0041] Characteristics such as disease, race, age, sex, or other biological or environmental factors that influence genetic traits can have a systematic impact on the distribution of raw gene or locus data:

[0042] (1) Disease status: including different types of malignant tumors, benign lesions, metabolic diseases, immune-related diseases and their different stages, subtypes or treatment response stages. Diseases can cause significant differences in gene or site signals between samples by altering cell composition, DNA methylation patterns or gene expression activity;

[0043] (2) Race or genetic background: refers to samples from different populations or genetic subgroups, whose differences in genotype and epigenetic modification may affect the level of site methylation, the frequency of single nucleotide variants (SNPs) and the distribution characteristics of detection signals;

[0044] (3) Age stage: The age difference of the individuals from which the sample is obtained can change the site signal through mechanisms such as tissue development status, degree of cell senescence and methylation drift, especially in the regions of development-related or aging-related genes;

[0045] (4) Sex characteristics: including differences in sex chromosome composition (such as XX and XY), sex hormone levels and related signaling pathway activity, which can affect X chromosome inactivation, sex-biased methylation patterns, etc.

[0046] (5) Tissue origin: refers to the type of biological material from which the sample is taken, such as blood, saliva, urine, tumor tissue, or peripheral tissue. Different tissues have specific cellular compositions and chromosome open states, which directly affect the distribution of site signal intensity and background noise;

[0047] (6) Lifestyle: including long-term behavioral factors such as smoking, drinking, diet, and exercise habits, which can indirectly change methylation or expression patterns by affecting metabolic state, oxidative stress level or inflammatory response pathways;

[0048] (7) Environmental exposure factors: including air pollution, heavy metals, drug exposure, radiation levels and occupational chemical exposure, which can cause systematic shifts in gene / locus signals through DNA damage repair pathways, epigenetic reprogramming or immune stress responses.

[0049] The factors mentioned above can act individually or in combination, leading to differences in the distribution, amplitude, and noise structure of the original signals of genes or loci among samples, which in turn affects the standardization of input data, feature selection, and classification accuracy during model construction.

[0050] The method described in this invention includes screening of differential features at multiple omics levels, specifically comprising:

[0051] (1) Transcriptome level: screening of differentially expressed genes, screening of coding genes for differentially expressed proteins, screening of differentially expressed miRNAs, and screening of non-coding RNAs based on differential expression abundance;

[0052] (2) Methylation and other epigenetic aspects: screening for differential methylation sites, screening for differential hydroxymethylation modification sites, screening for differential acetylation modification sites, screening for genes encoding proteins with differential phosphorylation ratios, and screening for genes with differential epigenetic modification intensity;

[0053] (3) At the genomic level: screening for differentially mutated genes, screening for differentially mutated sites, screening for differentially mutated copy number variant genes, screening for differentially mutated copy number structural variant genes, and screening for differentially frequent single nucleotide polymorphism (SNP) sites;

[0054] (4) Metabolic and pathway level: screening of genes related to metabolites with different concentrations, screening of genes encoding metabolic pathways with different activity, and screening of genes encoding metabolites with different contribution.

[0055] (5) Immune and microenvironment level: screening of coding genes related to differential infiltration ratio of immune cells, screening of coding genes related to differential concentration of cytokines, and screening of coding genes related to differential abundance of microbial communities;

[0056] (6) Molecular interaction and regulation level: screening of genes encoding proteins with different interaction strengths and screening of corresponding genes or sites for differential inhibition efficiency targets.

[0057] The sample types selected for gene / locus feature screening are selected from tumors, genetic abnormalities, metabolic abnormalities, diseases caused by drug-resistant pathogens, autoimmune diseases, neurodegenerative diseases, cardiovascular diseases, infectious diseases, inflammatory diseases, rare diseases, endocrine diseases, hematological diseases, developmental disorders, mental illnesses, skin diseases, degenerative musculoskeletal diseases, or diseases caused by environmental exposure, or subcategories under the above classifications.

[0058] The method described in this invention is characterized in that the bioinformatics database is selected from one or more of the following: gene expression and transcriptome database, genetic variation and population frequency database, cancer genome database, genotype-tissue expression database, non-coding RNA database, metabolic pathway and signaling pathway database, protein and interaction database, and cancer somatic mutation database.

[0059] The bioinformatics database used in the method described in this invention is selected from one or more of GEO, ENA, SRA, dbSNP, TCGA, gnomAD, GTEx, ICGC, PRIDE, miRBase, NONCODE, TarBase, miRGator, KEGG, STRING, and COSMIC.

[0060] GEO: A gene expression database maintained by NCBI, storing microarray and sequencing data submitted by research institutions worldwide, covering various research fields including oncology and non-oncology. Data formats include GEO Datasets (by experiment) and GEO Profiles (by gene).

[0061] ENA: A sequence archive operated by the European Bioinformatics Institute (EBI), similar to NCBI's SRA, storing raw high-throughput sequencing data (such as RNA-seq, WGS).

[0062] SRA: NCBI's Sequence Read Archive, stores raw sequencing data (such as FASTQ files) and supports the public sharing of high-throughput sequencing data.

[0063] dbSNP: NCBI's Single Nucleotide Polymorphism (SNP) database, which contains genetic variation information of humans and other species, and is often used in GWAS research.

[0064] gnomAD: A genome-wide / exome variation frequency database that integrates data from large-scale sequencing projects to provide allele frequencies in healthy individuals for assessing the pathogenicity of rare variants.

[0065] TCGA: Cancer Genome Atlas, contains multi-omics data (mutations, expression, methylation, etc.) and clinical information on 33 types of tumors, and is a core resource for cancer research.

[0066] ICGC: International Cancer Genome Consortium, covering 50 types of cancer, integrating TCGA data and providing broader global cohort information.

[0067] GTEx: A database that studies the relationship between genotype and tissue expression. It analyzes eQTLs (quantitative trait loci) from autopsy samples to reveal the impact of genetic variations on gene expression.

[0068] PRIDE: A proteomics data archive that stores raw data and analysis results from mass spectrometry experiments, supporting protein identification and quantification studies.

[0069] miRBase: a miRNA sequence and annotation database that provides information on mature miRNAs, precursor sequences, and target genes.

[0070] NONCODE: A comprehensive database of non-coding RNAs (lncRNA, circRNA, etc.) that includes functional annotations and expression profiles.

[0071] TarBase is a database for experimental validation of miRNA target genes, containing literature-supported evidence of miRNA-mRNA interactions.

[0072] miRGator: a toolkit for predicting miRNA function, providing target gene prediction, pathway analysis, and disease association information.

[0073] KEGG: A database of metabolic and signaling pathways.

[0074] STRING: A protein interaction prediction platform.

[0075] COSMIC: Cancer Somatic Mutation Database.

[0076] Another objective of this invention is to provide a highly specific gene / locus combination screening algorithm. Its core consists of a heuristic search algorithm and a scoring system to assess the likelihood that a target gene / locus combination exceeds the current optimal gene / locus combination. The implementation of the functions described in this invention relies on computers, computer programs such as Python, and operators with corresponding computer programming and data processing capabilities.

[0077] Another object of the present invention is to provide the application in the screening of gene / locus combinations based on heuristic search, including differential expression gene screening, differential methylation gene screening, differential mutation gene screening, differential copy number variation gene screening, coding gene screening of differentially expressed proteins, differential expression miRNA screening, differential methylation site screening, differential mutation site screening, differential hydroxymethylation modification site screening, differential acetylation modification site screening, single nucleotide polymorphism site screening with differential frequencies, structural variation gene screening with differential copy numbers, differential epigenetic modification intensity gene screening, non-coding RNA screening with differential expression abundances, coding gene screening of differentially phosphorylated ratio proteins, coding gene screening of differentially interacting intensity proteins, related gene screening of differentially concentrated metabolites, differential active metabolic pathway gene screening, coding gene screening of differentially infiltrated ratio immune cells, coding gene screening of differentially concentrated cytokines, coding gene screening of differentially abundant flora, coding gene screening of differentially contributing metabolites, corresponding gene / locus screening of differentially inhibited efficiency targets.

[0078] For the method of the present invention, preferably, the discretization thresholds are x1, x2, x3, x4, where 0 < x1 < x2 < x3 < x4 < 1, and the discretization thresholds x1, x2, x3, x4 divide the closed real number interval from 0 to 1 into discontinuous regions [0, x1], [x2, x3], [x4, 1].

[0079] Preferably, in step (4), the discretization thresholds are adjusted, and one or more of the values of x1...x w are increased or decreased by 0.01 to 0.05. More preferably, it is 0.02.

[0080] Preferably, the discretization thresholds are adjusted to magnify each interval.

[0081] Preferably, the following updates are respectively made to the four discretization thresholds: +0.02, -0.02, +0.02, -0.02.

[0082] A specific example of the present invention has the following parameters:

[0083] Discretization thresholds: The discretization thresholds are initially defined as 0.2, 0.4, 0.5, and 0.7. The discretization thresholds divide the closed real number interval from 0 to 1 into three discontinuous regions [0, 0.2], [0.4, 0.5], [0.7, 1];

[0084] k: The number of CpG sites in the current combination, initially initialized to 2;

[0085] km: The CpG site number threshold, which is the maximum value of the number of genes / loci for target differentiating sample types, initially initialized to 5;

[0086] m: No improvement cycle count, which is the number of consecutive cycles when there is no improvement in the scoring process. The initial value is defined as 1.

[0087] ct: No-improvement loop threshold, which is the maximum number of consecutive loops without any improvement during the scoring process, initialized to 50000;

[0088] n: The number of CpG sites to replace. This is the number of genes in the optimal combination of CpG sites that need to be replaced after reaching the no-improvement cycle threshold. The initial value is defined as 1.

[0089] Another objective of this invention is to provide a highly specific combination of genes / sites obtained by screening using the method described in this invention.

[0090] Another objective of this invention is to provide the application of the highly specific gene / site combination described herein in the preparation of sample type identification reagents.

[0091] The sample type identification reagents are disease diagnostic reagents, pathological testing reagents, pathogen detection reagents, tissue or cell type identification reagents, molecular marker detection reagents, biological sample identification reagents, or other testing reagents used in clinical practice, scientific research, and precision medicine.

[0092] Advantages of this invention:

[0093] This invention establishes a computational gene / locus combination screening method for predicting sample types that is highly specific, highly sensitive, and low-cost. The constructed heuristic search-based gene / locus screening algorithm simplifies the gene / locus combination screening process while maintaining the rationality and correctness of the optimization process. This ensures the model optimization process has the correct direction, making the final screening process more specific and efficient, significantly increasing the quality of sample type differentiation, and reducing the number of genes / locus in the gene / locus combination. The algorithm logic described in this invention forces the optimization of the gene / locus screening process towards minimizing the number of genes / locus and maximizing the differentiation of sample types, thereby making the screened gene / locus combinations more efficient and specific. Attached Figure Description

[0094] To more clearly illustrate the technical solution of the present invention, the accompanying drawings used in the embodiments of the present invention will be briefly described below.

[0095] Obviously, the accompanying drawings described below are only some of the drawings of the embodiments of the present invention. For those skilled in the art, other drawings can be obtained from these drawings without creative effort, but these other drawings are also within the scope of the drawings required for the embodiments of the present invention.

[0096] Figure 1This is a flowchart of the gene / site combination screening method of the present invention. Detailed Implementation

[0097] To make the objectives, technical solutions, beneficial effects, and significant advancements of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings provided in the embodiments of the present invention.

[0098] Obviously, all the embodiments described are only some embodiments of the present invention, and not all embodiments; based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0099] Those skilled in the art can understand the specific meaning of the above terms in this invention according to the specific circumstances.

[0100] It should also be noted that the following specific embodiments can be combined with each other, and the same or similar concepts or processes may not be described again in some embodiments.

[0101] The technical solution of the present invention will now be described in detail with reference to specific embodiments.

[0102] Example 1: An example application of the site selection process described in this invention.

[0103] 1. Data Preparation

[0104] The purpose of this embodiment is to demonstrate that the gene / site screening algorithm based on heuristic search described in this invention can screen out better site combinations. Specifically, this embodiment uses partial methylation sequencing data from a real disease-related database (TCGA database) as raw data and performs corresponding data processing. Then, the gene / site screening algorithm based on heuristic search described in this invention is used to find the CpG site combination with the best cancer (sample type) specificity, and the ability of the screened CpG site combination to distinguish different cancers is verified by a machine learning model.

[0105] Specifically, methylation profiles of genomic DNA from bladder urothelial carcinoma (TCGA-BLCA), colorectal cancer (TCGA-COAD & TCGA-READ), invasive breast cancer (TCGA-BRCA), esophageal cancer (TCGA-ESCA), and normal samples were obtained from the TCGA database. All selected data were sequencing data (IDAT format) from an Illumina 450k methylation chip, which contains methylation detection of 485,577 CpG sites. In this embodiment, a total of 412 bladder urothelial carcinoma samples, 393 colorectal cancer samples, 793 invasive breast cancer samples, 185 esophageal cancer samples, and 371 normal samples were collected. After quality control screening using the R package, 373,289 CpG sites were selected for analysis in this embodiment. This embodiment uses the beta value (methylation rate) after bioinformatics analysis (R language) as CpG site data, with the aim of finding s CpG sites that can effectively distinguish sample types (cancer). Specifically, the goal is to find the CpG site combination that achieves the greatest sample type discrimination ability while minimizing s.

[0106] 2. Screening for the optimal CpG site combination

[0107] The experiment was conducted in two groups: one group used an exhaustive search method to select the optimal CpG site combination, and the other group used the algorithm described in this invention to select the optimal CpG site combination with the minimum number of sites. In order to enable the exhaustive search method to analyze with a limited number of CpG sites, the number of CpG sites in the final target CpG site combination was limited to between 2 and 5.

[0108] Method 1: Exhaustive search

[0109] The value of 's' is set sequentially from 2 to 5. Then, 's' combinations of CpG sites are selected from all CpG sites. A neural network model is used to analyze the ability of each CpG site combination to distinguish different cancer types. The ability of each combination to distinguish all test set samples is evaluated sequentially. Finally, the optimal CpG site combination is selected.

[0110] For the exhaustive method, after each random selection of CpG site combinations, the methylation rate data of all independent samples corresponding to all sample types for all randomly selected CpG sites are split. The independent samples for each sample type are divided into training, validation, and test datasets, with the training dataset accounting for 70%, the validation dataset for 10%, and the test dataset for 20%. A neural network model is trained using the training and validation datasets. This neural network model contains an input layer, a hidden layer, and an output layer. The input layer contains 2-5 neurons, each receiving methylation rate data from 2-5 CpG sites. The hidden layer is a fully connected layer with 8 neurons and the ReLU activation function. The output layer contains 5 neurons, each corresponding to the probability of a different sample type, and the Softmax activation function is used. The output layer with the largest value is selected as the final predicted sample type result. The learning rate of the neural network model was set to 0.0001, and the loss function was Sparse Categorical Crossentropy. L2 regularization was used in each fully connected layer to avoid overfitting, with the L2 parameter set to 0.0001. Early stopping was implemented using Keras callbacks, with the early stopping condition being that the loss function value no longer decreased after five consecutive iterations. The training epochs were set to 100, and the number of input samples for each parameter update was set to 20. The labels during training corresponded to specific sample types. After training and building the model, it was used to predict the sample types of the test dataset. The model performance was assessed by comparing the predicted types with the actual sample types. The performance of each CpG site combination was exhaustively compared, and the CpG site combination with the best performance was selected as the final CpG site combination output using the exhaustive method.

[0111] In actual program execution, the exhaustive method requires performing sample type discrimination analysis on the 69672152116 combinations when s = 2, and on the 8.65 × 10⁻⁶ combinations when s = 3. 15 The ability of differentiating sample types was analyzed using various combinations, and the result was 8.07 × 10⁻⁶ when s = 4. 20 The ability of differentiating sample types was analyzed using various combinations, and the result was 6.02 × 10⁻⁶ when s = 5. 25 The ability of various combinations to distinguish sample types is analyzed. Even if the construction and prediction of a single neural network model are at the microsecond level, the exhaustive search method would take hundreds of millions of years to complete. Even with parallel computing using high-performance computers, it would be impossible to screen for the optimal combination of CpG sites.

[0112] Method 2: The CpG site screening method described in this invention

[0113] 1. Site screening, including the following screening steps:

[0114] (1) Initialize and define hyperparameters

[0115] Discretization threshold: Initially defined as 0.2, 0.4, 0.5, and 0.7, the discretization threshold divides the closed real number interval from 0 to 1 into three discontinuous regions: [0, 0.2], [0.4, 0.5], and [0.7, 1].

[0116] k: The number of CpG sites in the current combination, initialized to 2;

[0117] km: CpG site number threshold, the maximum number of genes / sites that can distinguish sample types, initialized to 5;

[0118] m: No improvement cycle count, which is the number of consecutive cycles when there is no improvement in the scoring process. The initial value is defined as 1.

[0119] ct: No-improvement loop threshold, which is the maximum number of consecutive loops without any improvement during the scoring process, initialized to 50000;

[0120] n: The number of CpG sites to be replaced, representing the number of CpG sites in the optimal CpG site combination to be replaced. The initial value is defined as 1.

[0121] (2) Discretization

[0122] For methylation microarray analysis data, the methylation rate is already within the range of [0, 1], so no further normalization is needed. Only the methylation data of all independent samples corresponding to each CpG site for each cancer type needs to be averaged. The data is divided according to the initially defined discretization thresholds. Data in the range [0, 0.2] is assigned the label a1, data in the range [0.4, 0.5] is assigned the label a2, and data in the range [0.7, 1] is assigned the label a3. Data not falling within any of the above discretization threshold ranges is assigned an empty label. Each CpG site corresponds to one label for each cancer type. CpG sites with identical labels for all sample types and all CpG sites containing empty labels are removed to obtain a standard CpG site dataset.

[0123] (3) Scoring is based on a heuristic search loop algorithm, specifically including:

[0124] (i) Randomly select k CpG sites from the CpG site standard dataset as the initial best CpG site combination, and obtain the label spectra formed by each label for each sample type corresponding to the best CpG site combination. Compare the label spectra of different sample types, and calculate the number of sample types with non-repeating label spectra, which is defined as the best score;

[0125] (ii)Calculate the number of sample types with non-repeating label spectra for the remaining k - 1 CpG site combinations after removing any one CpG site from the current best CpG site combination, which is defined as the CpG site score of the removed CpG site. A total of k CpG site scores are obtained. Sort the CpG site scores from largest to smallest, select the top n CpG sites with the highest rankings to form the CpG site index to be replaced. Re-select n CpG sites from the standard dataset after removing the CpG sites from the best CpG site combination to replace the CpG site index in the current best CpG site combination, construct a candidate CpG site combination, and obtain the label spectra formed by each label for each sample type corresponding to the candidate CpG site combination. Compare the label spectra of different sample types, and calculate the number of sample types with non-repeating label spectra, which is defined as the candidate score;

[0126] (iii)Scoring loop:

[0127] If the candidate score is greater than the best score, replace the current candidate CpG site combination with the best CpG site combination, replace the candidate score with the best score, refer to step (ii) to construct a new candidate CpG site combination and calculate a new candidate score, and repeat the above operations until the updated candidate score is less than or equal to the best score replaced in the previous round;

[0128] If the candidate score is less than or equal to the best score, keep the best CpG site combination unchanged, keep the best score unchanged, keep the CpG site scores unchanged, keep the CpG site index to be replaced unchanged, and refer to step (ii) to construct a new candidate CpG site combination and calculate a new candidate score, repeat the above operations, and update the m value by adding 1 successively with each increase in the number of loops until the new candidate score is greater than the best score;

[0129] (4)When the scoring loop in step (iii) reaches that the updated m value is greater than the ct value, perform parameter update:

[0130] (a)If n < k - 1, save the current best CpG site combination and the corresponding hyperparameters, reset the non-improved loop number m to 1, update the n value by adding 1, keep the best CpG site combination unchanged, keep the best score unchanged, keep the gene site scores unchanged, update the CpG site index to be replaced according to the updated n value, refer to step (ii) to construct a new candidate CpG site combination and calculate a new candidate score, and enter the scoring loop (iii);

[0131] (b) If n = k - 1, adjust the discretization thresholds and update the four discretization thresholds as follows: +0.02, -0.02, +0.02, -0.02, that is, enlarge each methylation rate interval to cover more available CpG sites. Then reset n to the initial value of 1, and reset the no-improvement cycle number m to 1.

[0132] Repeat steps (2) and (3) according to the new discretization threshold and update the parameters using steps (a) and (b) until the discretization threshold cannot be adjusted further. That is, when the four discretization thresholds are 0.3, 0.3, 0.6, and 0.6 respectively, if the current k is not greater than km, update the k value by 1, reset the discretization threshold to the initial discretization threshold, reset n to 1, reset m to 1, and repeat steps (2) and (3) and update the parameters using steps (a) and (b).

[0133] (5) Loop termination condition

[0134] When all loops and parameters are updated to the point where k is greater than km, exit all loops.

[0135] (6) Compare each of the best CpG site combinations and their corresponding hyperparameters saved in the loop, and select the CpG site combination with the largest best score as the final CpG site combination. If multiple best scores are the same, calculate the distance between each adjacent region of the closed real number interval from 0 to 1 divided by the current discretization threshold of each candidate CpG site combination and calculate the mean. Select the CpG site combination with the largest mean and / or the fewest CpG sites as the final CpG site combination.

[0136] 2. Assessment of the ability to distinguish between sample types

[0137] To demonstrate that the CpG sites screened by the method described in this invention have excellent cancer differentiation ability, a neural network model was used to analyze and verify the accuracy of the screened CpG sites in cancer differentiation.

[0138] For the method described in this invention, after completing the final selection of the optimal k CpG sites, the methylation rate data of all independent samples of all cancer types corresponding to the final selected k CpG sites are split. Each sample type's independent samples are divided into a training dataset, a validation dataset, and a test dataset, with the training dataset accounting for 70%, the validation dataset for 10%, and the test dataset for 20%. A neural network model is trained using the training and validation datasets. This neural network model includes an input layer, a hidden layer, and an output layer. The input layer contains k neurons, each receiving methylation rate data from the k CpG sites. The hidden layer is a fully connected layer with 8 neurons and the ReLU activation function. The output layer contains 5 neurons, each corresponding to the probability of a different sample type, and the Softmax activation function is used. The output layer with the largest value is selected as the final predicted sample type result. The learning rate of the neural network model was set to 0.0001, and the loss function was Sparse Categorical Crossentropy. L2 regularization was used in each fully connected layer to avoid overfitting, with the L2 parameter set to 0.0001. Early stopping was implemented using Keras callbacks, with the early stopping condition being that the loss function value no longer decreased after five consecutive iterations. The training epochs were set to 100, and the number of input samples for each parameter update was set to 20. The labels during training corresponded to specific sample types. After model construction, the model was used to predict the probability of each sample in the test dataset being classified as each sample type. The probabilities were sorted from highest to lowest, and the sample type with the highest probability was selected as the output. The prediction accuracy of this CpG site combination for the sample type (cancer) was calculated using the formula: the number of test set samples of that sample type (cancer type) accurately predicted as that sample type divided by the total number of test set samples of that sample type. The performance of CpG site screening based on this method was evaluated using prediction accuracy.

[0139] All neural network models were developed and run using Python. The high-level framework Keras was provided by TensorFlow for Python. The development environment consisted of Python 3.11.5 and TensorFlow 2.16.1, developed and run using the PyCharm Community Edition 2023.2.1 IDE. The platform for model development and operation was a computer running Windows 11 Pro 23H2, with a 13th Gen Intel(R) Core(TM) i7-13700KF 3.4 GHz processor and 64GB of dual-channel DDR4 Hynix memory.

[0140] 3. Experimental Results

[0141] The results of CpG site screening using the heuristic search-based gene / site screening algorithm described in this invention are shown in Table 1 below. Using the method described in this invention, a dual CpG site combination (chr3:37459206 and chr2:11484705) was screened within one hour. In the test dataset, this combination achieved a minimum accuracy of 94.59% for cancer identification.

[0142] Although exhaustive search is theoretically the only way to find the global optimal solution, it is limited by computational performance and complexity, making it impossible to perform a reasonable exhaustive analysis. Experimental results show that the cancer discrimination ability of the CpG site combination screened by the method described in this invention is excellent, that is, it can distinguish the most sample types with the fewest number of CpG sites, and at the same time, it has excellent performance. Therefore, the gene / site screening method described in this invention has strong efficiency and ability in screening the optimal CpG site combination.

[0143] Table 1: Statistical table of results for optimal CpG site selection using a heuristic search-based screening method.

[0144]

[0145] Example 2: An example application of the gene screening process described in this invention.

[0146] 1. Data Preparation

[0147] The purpose of this embodiment is to demonstrate that the gene / site screening algorithm based on heuristic search described in this invention can screen out better gene combinations. Specifically, this embodiment uses a portion of RNA sequencing data from a real disease-related database (TCGA database) as raw data and performs corresponding data processing. Then, the gene / site screening algorithm based on heuristic search described in this invention is used to find the gene / site combination with the best cancer (sample type) specificity, and the ability of the screened gene / site combination to distinguish different cancers is verified by a machine learning model.

[0148] Specifically, genomic RNA sequencing results of bladder urothelial carcinoma (TCGA-BLCA), colorectal cancer (TCGA-COAD & TCGA-READ), invasive breast cancer (TCGA-BRCA), esophageal cancer (TCGA-ESCA), and normal samples were obtained from the TCGA database. Illumina microarray sequencing was used, yielding RNA expression data for a total of 60,661 genes. All genes were analyzed using fpkm_uq_unstranded data. In this embodiment, 400 bladder urothelial carcinoma samples, 400 colorectal cancer samples, 800 invasive breast cancer samples, 185 esophageal cancer samples, and 400 normal samples were collected. The R package was used for data quality control and screening, ultimately determining the use of fpkm_uq_unstranded normalized data for 11,686 genes for analysis in this embodiment. The purpose of this embodiment is to find combinations of s genes that can effectively distinguish sample types. Specifically, the goal is to find gene combinations that achieve the greatest sample type discrimination ability while minimizing s.

[0149] 2. Screening for optimal gene combinations

[0150] The experiment was conducted in two groups: one group used an exhaustive search method to screen for the optimal gene combination, and the other group used the algorithm described in this invention to screen for the optimal gene combination with the minimum number of genes. In order to enable the exhaustive search method to analyze with a limited number of genes, the number of genes in the final target gene combination was limited to between 2 and 5.

[0151] Method 1: Exhaustive search

[0152] By sequentially setting s to values ​​from 2 to 5, we exhaustively select s gene combinations from all genes. A neural network model is then used to analyze the ability of each gene combination to distinguish different sample types, evaluating the ability of each combination to distinguish all test set samples. Finally, the optimal gene combination is selected through comparison.

[0153] For the exhaustive method, after each random selection of gene combinations, the gene expression level data of all independent samples corresponding to all cancer types obtained by the random selection are split. The independent samples of each sample type are divided into training dataset, validation dataset, and test dataset, with the training dataset accounting for 70%, the validation dataset accounting for 10%, and the test dataset accounting for 20%. A neural network model is trained using the training dataset and the validation dataset. The neural network model contains an input layer, a hidden layer, and an output layer. The input layer contains 2 to 5 neurons, which are used to receive 2 to 5 gene expression level data. The hidden layer is a fully connected layer with 8 neurons and the activation function is ReLU. The output layer contains 5 neurons, which correspond to the probabilities of different sample types and the activation function is Softmax. The output layer with the largest value is selected as the final predicted sample type result of the model. The learning rate of the neural network model was set to 0.0001, and the loss function was Sparse Categorical Crossentropy. L2 regularization was used in each fully connected layer to avoid overfitting, with the L2 parameter set to 0.0001. Early stopping was added using Keras callbacks, with the early stopping condition being that the loss function value no longer decreased after five consecutive iterations. The number of training epochs was set to 100, and the number of input samples for each parameter update was set to 20. The labels during training corresponded to specific sample types. After training and building the model, it was used to predict the sample types of the test dataset. The model performance was judged by comparing the predicted results with the actual sample types. The performance of each gene combination was exhaustively compared, and the gene combination with the best model performance was selected as the final gene combination output using the exhaustive method.

[0154] In the actual operation of the program, it is necessary to perform sample type discrimination ability analysis on 68,306,741 combinations when s = 2, and on 2.66 × 10⁻⁶ combinations when s = 3. 11 The ability of differentiating sample types was analyzed using various combinations, and the 7.77 × 10⁻⁶ samples were analyzed when s = 4. 14 The ability of differentiating sample types was analyzed using various combinations, and the result was 1.81 × 10⁻⁶ when s = 5. 18 The ability of various combinations to distinguish sample types is analyzed. Even if the construction and prediction of a single neural network model are at the microsecond level, the exhaustive search method would take tens of thousands of years to complete, and even using high-performance computers for parallel computing would not be able to complete the screening of the optimal gene combination.

[0155] Method 2: The gene combination screening method described in this invention

[0156] 1. Gene combination screening, including the following screening steps:

[0157] (1) Initialize and define hyperparameters

[0158] Discretization threshold: Initially defined as 0.2, 0.4, 0.5, and 0.7, the discretization threshold divides the closed real number interval from 0 to 1 into three discontinuous regions: [0, 0.2], [0.4, 0.5], and [0.7, 1].

[0159] k: The number of genes in the current combination, initialized to 2;

[0160] km: Gene count threshold, the maximum number of genes required to distinguish sample types, initialized to 5;

[0161] m: No improvement cycle count, which is the number of consecutive cycles when there is no improvement in the scoring process. The initial value is defined as 1.

[0162] ct: No-improvement loop threshold, which is the maximum number of consecutive loops without any improvement during the scoring process, initialized to 50000;

[0163] n: The number of genes to be replaced, the number of genes in the optimal gene combination to be replaced, initially defined as 1.

[0164] (2) Discretization

[0165] All raw data are normalized by performing a log transformation on the fpkm_uq_unstranded for each gene in each sample, using the operation log(X+1, 2) to obtain the transformed fpkm_uq_unstranded, where X is the original fpkm_uq_unstranded for any gene. Then, the log-transformed fpkm_uq_unstranded is normalized by dividing all fpkm_uq_unstranded by the maximum value of fpkm_uq_unstranded in that sample. Finally, the transformed fpkm_uq_unstranded of all independent samples for all cancer types are averaged. The processed data is divided according to the discretization threshold defined in the initialization. Each data point in the range [0, 0.2] is assigned a label a1, data in the range [0.4, 0.5] is assigned a label a2, and data in the range [0.7, 1] is assigned a label a3. Each gene has a label corresponding to each cancer type. Data points that do not fall within any of the above discretization threshold ranges are assigned an empty label. Genes that have the same label for all sample types and all genes containing empty labels are removed to obtain the standard gene dataset.

[0166] (3) Scoring is based on a heuristic search loop algorithm, specifically including:

[0167] (i) Randomly select k genes from the gene standard dataset as the initialized optimal gene combination, and obtain the label spectra formed by each label for each sample type corresponding to the optimal gene combination. Compare the label spectra of different sample types, and calculate the number of sample types with non-repeating label spectra, which is defined as the optimal score.

[0168] (ii) Calculate the number of sample types with non-repeating label spectra for the remaining k - 1 gene combinations after removing any one gene from the current optimal gene combination, which is defined as the gene score of the removed gene. A total of k gene scores are obtained. Sort the gene scores from largest to smallest, select the top n genes with the highest rankings to form the gene index to be replaced, reselect n genes from the standard dataset after removing the genes from the optimal gene combination to replace the gene index in the current optimal gene combination, construct a candidate gene combination, and obtain the label spectra formed by each label for each sample type corresponding to the candidate gene combination. Compare the label spectra of different sample types, and calculate the number of sample types with non-repeating label spectra, which is defined as the candidate score.

[0169] (iii) Scoring loop:

[0170] If the candidate score is greater than the optimal score, replace the current candidate gene combination with the optimal gene combination, replace the candidate score with the optimal score, construct a new candidate gene combination and calculate a new candidate score according to step (ii), and repeat the above operations until the updated candidate score is less than or equal to the optimal score replaced in the previous round.

[0171] If the candidate score is less than or equal to the optimal score, keep the optimal gene combination unchanged, keep the optimal score unchanged, keep the gene site scores unchanged, keep the gene index to be replaced unchanged, and construct a new candidate gene combination and calculate a new candidate score according to step (ii), repeat the above operations, and update the m value by adding 1 successively with each increase in the number of loops until the new candidate score is greater than the optimal score.

[0172] (4) When the m value updated in the scoring loop of step (iii) is greater than the ct value, perform parameter update:

[0173] (a) If n < k - 1, save the current optimal gene combination and the corresponding hyperparameters, reset the non-improved loop number m to 1, update the n value by adding 1, keep the optimal gene combination unchanged, keep the optimal score unchanged, keep the gene site scores unchanged, update the gene index to be replaced according to the updated n value, construct a new candidate gene combination and calculate a new candidate score according to step (ii), and enter the scoring loop (iii).

[0174] (b) If n = k - 1, adjust the discretization thresholds as follows: +0.02, -0.02, +0.02, -0.02, respectively. This enlarges the expression level range of each gene to cover more available genes. Reset n to its initial value of 1, and reset the no-improvement loop number m to 1.

[0175] Repeat steps (2) and (3) according to the new discretization threshold and update the parameters using steps (a) and (b) until the discretization threshold cannot be adjusted further. That is, when the four discretization thresholds are 0.3, 0.3, 0.6, and 0.6 respectively, if the current k is not greater than km, update the k value by 1, reset the discretization threshold to the initial discretization threshold, reset n to 1, reset m to 1, and repeat steps (2) and (3) and update the parameters using steps (a) and (b).

[0176] (5) Loop termination condition

[0177] When all loops and parameters are updated to the point where k is greater than km, exit all loops.

[0178] (6) Compare each of the best gene combinations and their corresponding hyperparameters saved in the loop, and select the genome combination with the largest best score as the final target gene combination. If multiple largest best scores are the same, calculate the distance between each adjacent region of the closed real number interval from 0 to 1 divided by the current discretization threshold of each candidate gene combination and calculate the mean. Select the genome combination with the largest mean and / or the fewest number of genes as the final target gene combination.

[0179] 2. Assessment of the ability to distinguish between sample types

[0180] To demonstrate that the gene combinations screened by the method described in this invention have excellent cancer differentiation ability, a neural network model was used to analyze and verify the accuracy of the screened gene combinations in differentiating cancers.

[0181] In the method described in this invention, after selecting the final optimal k genes, the gene expression data of all independent samples corresponding to all cancer types for the final selected k genes are split. Each independent sample type is divided into a training dataset, a validation dataset, and a test dataset, with the training dataset accounting for 70%, the validation dataset for 10%, and the test dataset for 20%. A neural network model is trained using the training and validation datasets. This neural network model includes an input layer, a hidden layer, and an output layer. The input layer contains k neurons, each receiving gene expression level data for the k genes. The hidden layer is a fully connected layer with 8 neurons and the ReLU activation function. The output layer contains 5 neurons, each corresponding to the probability of a different sample type, and the Softmax activation function is used. The output layer with the largest value is selected as the final predicted sample type result. The learning rate of the neural network model was set to 0.0001, and the loss function was Sparse Categorical Crossentropy. L2 regularization was used in each fully connected layer to avoid overfitting, with the L2 parameter set to 0.0001. Early stopping was implemented using Keras callbacks, with the early stopping condition being that the loss function value no longer decreased after five consecutive iterations. The training epochs were set to 100, and the number of input samples for each parameter update was set to 20. The labels during training corresponded to specific sample types. After model construction, the model was used to predict the probability of each sample in the test dataset being classified as each sample type. The probabilities were sorted from highest to lowest, and the sample type with the highest probability was selected as the output. The prediction accuracy of the gene combination for each sample type was calculated using the formula: the number of test dataset samples of that type accurately predicted as that type divided by the total number of test dataset samples of that type. The performance of gene screening based on this method was evaluated using the prediction accuracy.

[0182] All neural network models were developed and run using Python. The high-level framework Keras was provided by TensorFlow for Python. The development environment consisted of Python 3.11.5 and TensorFlow 2.16.1, developed and run using the PyCharm Community Edition 2023.2.1 IDE. The platform for model development and operation was a computer running Windows 11 Pro 23H2, with a 13th Gen Intel(R) Core(TM) i7-13700KF 3.4 GHz processor and 64GB of dual-channel DDR4 Hynix memory.

[0183] 3. Experimental Results

[0184] The results of gene screening using the heuristic search-based gene / site screening algorithm described in this invention are shown in Table 2 below. Using the method described in this invention, the dual gene combination (ENSG00000170469.11 and ENSG00000130985.17) was screened within one hour. In the test dataset, the lowest accuracy of this gene combination in identifying various cancers was as high as 91.89%.

[0185] Although exhaustive search is theoretically the only way to find the global optimal solution, it is limited by computational performance and complexity, making it impossible to perform a reasonable exhaustive analysis. Therefore, exhaustive search is not feasible in theory or engineering. Experimental results show that the gene combinations screened by the method described in this invention have excellent cancer discrimination ability, that is, they can distinguish the most sample types with the fewest number of genes, and at the same time, they have excellent performance. Therefore, the gene / site screening method described in this invention has strong efficiency and ability in screening the optimal gene combination.

[0186] Table 2: Statistical table of results for screening optimal gene combinations using a heuristic search-based screening method.

[0187]

Claims

1. A gene / locus screening method based on data discretization and heuristic search loop, characterized in that... Includes the following steps: (1) Initialize and define hyperparameters: Discretization threshold: Discretization threshold x1~x w 0 <x1<……<x w <1, Discretization threshold x1~x w Divide the closed real number interval from 0 to 1 into discontinuous regions [0, x1], ..., [x...]. w [1], w represents even numbers from 2 to 10; k: Number of genes / locus in the current combination; km: Gene / locus number threshold, which is the maximum number of genes / locus points that can be used to distinguish sample types; m: No improvement cycle count, which is the number of consecutive cycles when there is no improvement in the scoring process. The initial value is defined as 1. ct: No Improvement Cycle Threshold, which is the maximum number of consecutive cycles without any improvement during the scoring process; n: The number of genes / sites to be replaced. The number of genes / sites in the optimal combination of genes / sites to be replaced. The initial value is defined as 1. (2) Discretization process: Select several independent samples under each sample type from the bioinformatics database to obtain the original data of all genes / locus of all independent samples corresponding to each sample type. Normalize the original data that are not in the range [0, 1] so that all data are scaled to within 0 to 1. Then, take the average value of the data of all independent samples corresponding to each gene / locus under each sample type. The sample type is selected from the subclass of the same category of features of disease, race, age, sex or other biological or environmental factors that affect gene characteristics. The mean-normalized data is divided according to the initially defined discretization threshold. Each data point in the range [0, x1] is assigned a label a1, ..., and the data in the range [x... w The data within the range of 1] is assigned the label a (w / 2+1) Data that does not fall within any of the above discretization thresholds are assigned empty labels. Each gene / locus has one label for each sample type. Genes / locus that have the same label for all sample types and all genes / locus containing empty labels are removed to obtain a standard gene / locus dataset. (3) Scoring is based on a heuristic search loop algorithm, specifically including: (i) Randomly select k genes / locus from the standard gene / locus dataset as the optimal gene / locus combination for initialization, and obtain the label spectrum formed by each label for each sample type corresponding to the optimal gene / locus combination. Compare the label spectrum of different sample types, calculate the number of sample types with non-repeating label spectrum, and define it as the optimal score. (ii) Calculate the number of unique sample types in the label spectrum of the remaining k-1 gene / site combinations after removing any gene / site from the current best gene / site combination. This number is defined as the gene / site score of the removed gene / site. A total of k gene / site scores are obtained. The gene / site scores are sorted from largest to smallest. The top n genes with the highest rankings are selected to form the gene / site index to be replaced. From the standard dataset where the genes / sites in the best gene / site combination have been removed, n gene / sites are reselected to replace the gene / site index to be replaced in the current best gene / site combination. Candidate gene / site combinations are constructed, and the label spectrum of each sample type corresponding to the candidate gene / site combination is obtained. The label spectra of different sample types are compared, and the number of unique sample types in the label spectrum is calculated. This number is defined as the candidate score. (iii) Scoring Cycle: If the candidate score is greater than the best score, the current candidate gene / site combination is replaced with the best gene / site combination, and the candidate score is replaced with the best score. Refer to step (ii) to construct a new candidate gene / site combination and calculate a new candidate score. Repeat the above operation until the updated candidate score is less than or equal to the best score of the previous round of replacement. If the candidate score is less than or equal to the best score, keep the best gene / locus combination unchanged, keep the best score unchanged, keep the gene / locus score unchanged, keep the index of the gene / locus to be replaced unchanged, and refer to step (ii) to construct a new candidate gene / locus combination and calculate a new candidate score. Repeat the above operations, and update the value of m by adding 1 successively with each increase in the number of cycles until the new candidate score is greater than the best score; (4) When the score loop in step (iii) reaches the updated m value greater than the ct value, perform parameter update: (a) If n < k - 1, save the current best gene / locus combination and the corresponding hyperparameters, reset the non-improvement cycle number m to 1, update the value of n by adding 1, keep the best gene / locus combination unchanged, keep the best score unchanged, keep the gene locus score unchanged, update the index of the gene / locus to be replaced according to the updated value of n, refer to step (ii) to construct a new candidate gene / locus combination and calculate a new candidate score, and enter the score loop (iii); (b) If n=k-1, then adjust the discretization threshold by increasing and / or decreasing x1 to x w One or more values ​​in the set, and reset n to the initial value 1, and reset the no-improvement loop number m to 1; Repeat steps (2) and (3) according to the new discretization threshold and update the parameters using steps (a) and (b) until the discretization threshold cannot be further adjusted. That is, when the latest discretization threshold divides the closed interval from 0 to 1 into completely continuous sub-intervals, if the current k is not greater than km, update the value of k by adding 1, reset the discretization threshold to the initial discretization threshold, reset n to 1, reset m to 1, and repeat steps (2), (3) and update the parameters using steps (a) and (b); (5) Loop termination condition When all loops and parameter updates reach the current k greater than km, exit all loops; (6) Compare each best gene / locus combination saved in the loop and its corresponding hyperparameters, select the gene / locus combination with the largest best score value as the final target gene / locus combination. If there are multiple equal largest best score values, calculate the mean of the distances between adjacent regions of the closed real number interval from 0 to 1 divided by the current discretization threshold for each candidate gene / locus combination, and select the gene / locus combination with the largest mean and / or the fewest number of genes / loci as the final target gene / locus combination.

2. The method as described in claim 1, characterized in that... The original data is one or several of DNA sequencing data, RNA sequencing data, DNA methylation sequencing data, miRNA sequencing data, protein expression profiles, DNA copy number variation profiles, metabolomics data, genomics data, proteomics data, transcriptomics data, epigenomics data, immunomics data, pharmacomics data, systems biology data, phenomics data, lipidomics data, glycomics data, signalomics data, toxicomics data or evolutionary genomics data.

3. The method as described in claim 2, characterized in that... The raw data are selected from one or more of the following: methylation rate, gene expression level, mutation status, gene copy number, miRNA expression abundance, protein expression level, single nucleotide polymorphism frequency, structural variation copy number, epigenetic modification intensity, non-coding RNA expression abundance, phosphorylated protein ratio, protein-protein interaction intensity, metabolite concentration, metabolic pathway activity, immune cell infiltration ratio, cytokine concentration, microbial abundance, metabolite contribution, target inhibition efficiency, glucose metabolism concentration, or lipid metabolism concentration.

4. The method as described in claim 1, characterized in that, The other biological or environmental factors that influence gene characteristics are selected from genetic background, tissue origin, lifestyle, and environmental exposure factors.

5. The gene / site screening method based on heuristic search as described in claim 1, characterized in that, The screening includes differential feature screening at multiple omics levels, specifically including: (1) Transcriptome level: screening of differentially expressed genes, screening of coding genes for differentially expressed proteins, screening of differentially expressed miRNAs, and screening of non-coding RNAs based on differential expression abundance; (2) Methylation and other epigenetic aspects: screening for differential methylation sites, screening for differential hydroxymethylation modification sites, screening for differential acetylation modification sites, screening for genes encoding proteins with differential phosphorylation ratios, and screening for genes with differential epigenetic modification intensity; (3) At the genomic level: screening for differentially mutated genes, screening for differentially mutated sites, screening for differentially mutated copy number variant genes, screening for differentially mutated copy number structural variant genes, and screening for differentially frequent single nucleotide polymorphism (SNP) sites; (4) Metabolic and pathway level: screening of genes related to metabolites with different concentrations, screening of genes encoding metabolic pathways with different activity, and screening of genes encoding metabolites with different contribution. (5) Immune and microenvironment level: screening of coding genes related to differential infiltration ratio of immune cells, screening of coding genes related to differential concentration of cytokines, and screening of coding genes related to differential abundance of microbial communities; (6) Molecular interaction and regulation level: screening of genes encoding proteins with different interaction strengths and screening of corresponding genes or sites for differential inhibition efficiency targets.

6. The method as described in claim 5, characterized in that... The sample types for gene / locus feature screening are selected from subcategories under the categories of tumors, genetic abnormalities, metabolic disorders, diseases caused by drug-resistant pathogens, autoimmune diseases, neurodegenerative diseases, cardiovascular diseases, infectious diseases, inflammatory diseases, rare diseases, endocrine diseases, hematological diseases, developmental disorders, mental illnesses, skin diseases, degenerative musculoskeletal diseases, or diseases caused by environmental exposure.

7. The method as described in claim 1, characterized in that... The bioinformatics database is selected from one or more of the following: gene expression and transcriptome database, genetic variation and population frequency database, cancer genome database, genotype-tissue expression database, non-coding RNA database, metabolic pathway and signaling pathway database, protein and interaction database, and cancer somatic mutation database.

8. The method as described in claim 1, characterized in that... The bioinformatics database is selected from one or more of GEO, ENA, SRA, dbSNP, TCGA, gnomAD, GTEx, ICGC, PRIDE, miRBase, NONCODE, TarBase, miRGator, KEGG, STRING, COSMIC.

9. The method as described in claim 1, characterized in that... The initial gene / locus number threshold (km) is not less than log0. (Dis_C) Sam_C is an integer, where Dis_C is the number of discretization intervals, Sam_C is the number of sample types to be distinguished, and the no-improvement loop threshold ct is an integer greater than or equal to X / 5, where X is the number of all optional genes / sites.

10. The method as described in claim 1, characterized in that... The discretization thresholds are x1, x2, x3, x4, where 0 < x1 < x2 < x3 < x4 < 1, and the discretization thresholds x1, x2, x3, x4 divide the closed real number interval from 0 to 1 into discontinuous regions [0, x1], [x2, x3], [x4, 1].

11. The method as described in claim 1, characterized in that... Step (4) Adjust the discretization threshold, setting x1...x w One or more values ​​in the range may increase or decrease by 0.01 to 0.

05.

12. A highly specific gene / site combination, characterized in that... It is obtained by screening using the method according to any one of claims 1 to 11.

13. Use of the highly specific gene / locus combination according to claim 12 in the preparation of a sample type identification reagent.

14. The application as described in claim 13, characterized in that... The sample type identification reagent is a disease diagnosis reagent, a pathological detection reagent, a pathogen detection reagent, a tissue or cell type identification reagent, a molecular marker detection reagent, a biological sample identification reagent, or other detection reagents related to clinical, scientific research and precision medicine.

Images