WO2013026411A1 - 单细胞分类方法、基因筛选方法及其装置 - Google Patents

单细胞分类方法、基因筛选方法及其装置 Download PDF

Info

Publication number
WO2013026411A1
WO2013026411A1 PCT/CN2012/080578 CN2012080578W WO2013026411A1 WO 2013026411 A1 WO2013026411 A1 WO 2013026411A1 CN 2012080578 W CN2012080578 W CN 2012080578W WO 2013026411 A1 WO2013026411 A1 WO 2013026411A1
Authority
WO
WIPO (PCT)
Prior art keywords
single cell
genotype
snp
reads
cell
Prior art date
Application number
PCT/CN2012/080578
Other languages
English (en)
French (fr)
Inventor
徐讯
鲍莉
何伟明
侯勇
陶晔
Original Assignee
深圳华大基因科技有限公司
深圳华大基因研究院
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Family has litigation
First worldwide family litigation filed litigation Critical https://patents.darts-ip.com/?family=47745956&utm_source=google_patent&utm_medium=platform_link&utm_campaign=public_patent_search&patent=WO2013026411(A1) "Global patent litigation dataset” by Darts-ip is licensed under a Creative Commons Attribution 4.0 International License.
Application filed by 深圳华大基因科技有限公司, 深圳华大基因研究院 filed Critical 深圳华大基因科技有限公司
Priority to US14/239,650 priority Critical patent/US20140206006A1/en
Priority to EP12825673.2A priority patent/EP2749655B2/en
Publication of WO2013026411A1 publication Critical patent/WO2013026411A1/zh

Links

Classifications

    • C—CHEMISTRY; METALLURGY
    • C12—BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
    • C12Q—MEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
    • C12Q1/00—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions
    • C12Q1/68—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions involving nucleic acids
    • C12Q1/6813—Hybridisation assays
    • C12Q1/6827—Hybridisation assays for detection of mutation or polymorphism
    • C—CHEMISTRY; METALLURGY
    • C12—BIOCHEMISTRY; BEER; SPIRITS; WINE; VINEGAR; MICROBIOLOGY; ENZYMOLOGY; MUTATION OR GENETIC ENGINEERING
    • C12Q—MEASURING OR TESTING PROCESSES INVOLVING ENZYMES, NUCLEIC ACIDS OR MICROORGANISMS; COMPOSITIONS OR TEST PAPERS THEREFOR; PROCESSES OF PREPARING SUCH COMPOSITIONS; CONDITION-RESPONSIVE CONTROL IN MICROBIOLOGICAL OR ENZYMOLOGICAL PROCESSES
    • C12Q1/00—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions
    • C12Q1/68—Measuring or testing processes involving enzymes, nucleic acids or microorganisms; Compositions therefor; Processes of preparing such compositions involving nucleic acids
    • C12Q1/6869—Methods for sequencing
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B20/00—ICT specially adapted for functional genomics or proteomics, e.g. genotype-phenotype associations
    • G16B20/20—Allele or variant detection, e.g. single nucleotide polymorphism [SNP] detection
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16B—BIOINFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR GENETIC OR PROTEIN-RELATED DATA PROCESSING IN COMPUTATIONAL MOLECULAR BIOLOGY
    • G16B30/00—ICT specially adapted for sequence analysis involving nucleotides or amino acids
    • G16B30/10—Sequence alignment; Homology search

Definitions

  • the present invention relates to bioinformatics, and more particularly to single cell classification and gene screening methods and apparatus for use in the methods. Background technique
  • Single cell classification can be effectively applied to the study of various stem cell differentiation processes, such as the differentiation of tumor stem cells, embryonic stem cells, and hematopoietic stem cells. It is necessary to screen stem cells at different stages of differentiation and perform various stem cell tests. In the study of drug resistance, it is necessary to accurately classify cells in different periods of administration, thereby further analyzing the drug resistance and drug resistance genes of the cell subset, for example, multidrug resistance and multidrug resistance in cancer patients. Study of the relationship between drug genes and drug abuse, drug tolerance, and drug dependence. Similarly, in the screening of drug target genes, since the drug interacts with cells, especially sensitive cells, it will cause a series of changes in the external morphology of the cells and the internal normal metabolic processes. Therefore, screening sensitive cells is the key first step.
  • the present invention aims to solve at least one of the technical problems existing in the prior art.
  • probability file of genotype refers to a set of values of the posterior probabilities of possible genotypes in the target region of the sample calculated using the Bayesian formula using SNP detection software; when using SNP detection When the software is SOAPsrp, the "generic possibility file” generated is the CNS file.
  • genotype file refers to the selection of the genotype with the highest probability among the above “generic type probability files” as the consistent genotype of each cell, based on the reference genome SNP dataset position information, extracting each The corresponding site of the cell genotype, the obtained genotype set of the population SNP at the corresponding site of each cell.
  • the invention proposes a single cell classification method.
  • the single cell classification method of the present invention comprises: sequencing a whole genome of a plurality of single cell samples from the same population separately to obtain reads from each single cell sample; respectively, from each single The reads of the cell samples are compared with the reference genome sequence, and the reads are subjected to data filtering; based on the filtered reads, the uniform genotype of each single cell sample is determined, wherein the uniform genotypes of all single cell samples are a set of SNP data for the population; for each of the single cells, based on the SNP data set of the population, determining a genotype corresponding to the position of the reference genomic SNP data set, each cell; and selecting A SNP site associated with a cell mutation, and the single cell is classified based on the genotype of the single cell at that site.
  • a single-generation genome can be analyzed and studied by a bioinformatics method using a new generation sequencing technology (NGS), and a subpopulation of cells (or microparticles) can be collected for further follow-up research.
  • NGS new generation sequencing technology
  • avoiding the labeling operation of cells effectively solves the problem that some cell subpopulations in the traditional single cell classification method do not have corresponding specific markers and cannot be classified; on the other hand, comprehensively and completely analyze single cells.
  • the genetic variation information of the genome, and multiple control samples can be set, which greatly improves the accuracy of cell subgroup classification.
  • the above single cell classification method may further have the following additional technical features:
  • the sequencing is performed using a second-generation or third-generation sequencing platform.
  • the data filtering criteria are: When there are multiple pairs of repeated double-end reads, the sequences of the pairs of reads are completely identical. A pair of reads are randomly selected, and other repeated double-end reads in the plurality of pairs of reads are removed; and/or reads that are not uniquely aligned to the reference genome sequence are removed.
  • determining the uniform genotype of each single cell based on the filtered reads further comprises: determining a likelihood of each single cell sample being genotyped in the target region based on the filtered reads; Based on the likelihood that all single cell samples will be genotyped in the target region, a pseudogenome containing each site for all samples is determined; and the most probable genotype is selected from the pseudogenome as the consensus genotype for each single cell sample.
  • selecting a SNP site associated with a cell mutation further removes at least one of the following from the population SNP data set: a non-community SNP site, a site in which a heterozygous loss exists, and SNP sites have been published.
  • the whole genome of at least one of the plurality of single cell samples is prior to sequencing Whole genome amplification treatment is performed, wherein removing the presence of heterozygous lost sites further includes removing sites that meet the following conditions: for samples that have not been whole genome amplified, the sequencing results are heterozygous sites; In the genomic amplified sample, the number of samples at the same site, heterozygous site and data loss is greater than or equal to the number of samples amplified by whole genome minus 3.
  • the genotype corresponding to each cell further The screening of the SNP data set is performed according to the following criteria: the quality value of the consistent genotype of each locus is not less than 20, the rank test p value is not less than 1%; and the SNPs for heterozygous variation: the best base
  • the quality of the sequencing type is not less than 20, the sequencing depth is not less than 6, the sub-base type sequencing quality value is not less than 20, the sequencing depth is not less than 2, and the ratio of the two genotype sequencing depths is in the range of 0.2-5.
  • the method further comprises the steps of: extracting each cell sample information and rejecting the objectioned cells.
  • after classifying the single cells further comprising: determining the classified population based on the classification result, and calculating a statistic of all SNP sites of each gene in each type of population, optionally The obtained statistic is subjected to a difference test to obtain a test value; a gene having the highest statistic or test value is selected as a gene associated with the cell mutation.
  • the invention provides a single cell sorting device.
  • the single cell classification device comprises: a data filtering module, the data filtering module being adapted to compare the reads from each single cell sample with a reference genome sequence, respectively, and perform the reading on the reads Data filtering, wherein the reading of each single cell sample is obtained by separately sequencing a whole genome of a plurality of single cell samples; a genotype determining module, wherein the genotype determining module is adapted to be based on filtered reads, Determining a consistent genotype for each single cell sample, wherein a uniform genotype of all single cell samples constitutes a SNP data set for said population; a genotype file extraction module, said genotype file extraction module being adapted for each of said Single cells, based on the SNP data set of the population, determining a genotype corresponding to each of the cells corresponding to the position of the reference genome SNP data set; a classification module, the classification module being adapted to be based on
  • the aforementioned single cell classification method can be effectively implemented by the apparatus.
  • a single-generation genome can be analyzed and studied by a bioinformatics method using a next-generation sequencing technology (NGS), and a subpopulation of cells (or microparticles) can be collected for further follow-up research.
  • NGS next-generation sequencing technology
  • avoiding the labeling operation of cells effectively solves the problem that some cell subpopulations in the traditional single cell classification method do not have corresponding specific markers and cannot be classified; on the other hand, comprehensively and completely analyze single cells.
  • the genetic variation information of the genome, and multiple control samples can be set, which greatly improves the accuracy of cell subgroup classification.
  • the single cell sorting device may also have the following additional technical features:
  • the data filtering module is adapted to perform data filtering based on the following criteria: when there are multiple pairs of repeated double-end reads, randomly selecting a pair of reads when the sequences of the pairs of reads are completely identical, and Other repeating double-end reads in the plurality of pairs of reads; and/or removing reads from the reference genome sequence not uniquely aligned
  • the genotyping module is adapted to determine a uniform genotype of each of the single cells by: determining, based on the filtered reads, a genotype of each single cell sample in a target region Possibility; determining the pseudogenome containing each locus of all samples based on the likelihood of all single-cell samples in the genotype of the target region; and selecting the most probable genotype from the pseudogenome as each single-cell sample Consistent genotype.
  • the classification module is adapted to remove at least one of the following from the population SNP data set to select a SNP site associated with a cell mutation: a non-community SNP site, a heterozygous loss Sites, and published SNP sites.
  • the whole genome of at least one of the plurality of single cell samples is subjected to whole genome amplification processing prior to performing sequencing, wherein the classification module is adapted to remove a site that meets the following conditions, In order to remove sites with loss of heterozygosity: for samples that have not been amplified by whole genome, the sequencing results are heterozygous sites; and in the genome-wide amplified samples, at the same site, heterozygous sites and data The number of samples lost is greater than or equal to the number of samples that have been amplified by the whole genome minus 3.
  • the genotype file extraction module is adapted to screen the SNP data set according to the following criteria: the quality value of the consistent genotype of each locus is not less than 20, and the rank test p value Not less than 1%; and SNPs for heterozygous variation: the best base type sequencing quality value is not less than 20, the sequencing depth is not less than 6, the sub-base type sequencing quality value is not less than 20, the sequencing depth is not less than 2, and The ratio of the depth of sequencing of the two genotypes is in the range of 0.2-5.
  • the classification module is further adapted to extract each cell sample information and reject the objectioned cells.
  • the screening module is adapted to: determine the classified population based on the classification result, and calculate a statistic of all SNP sites of each gene in each type of population, optionally The obtained statistic is subjected to a difference test to obtain a test value; and a gene having the highest statistic or test value is selected as a gene associated with the cell mutation.
  • the invention provides a method of genetic screening.
  • the method comprises the steps of: classifying cells to obtain a classified subpopulation, and calculating a statistic for all SNP loci for each gene in each subpopulation, optionally for the obtained The statistic is used for the difference test to obtain the test value; the gene with the highest statistic or test value is selected as the gene associated with the cell mutation.
  • pre-classifying cells for example, according to predetermined criteria, such as cancer cells and cancer cells, or other cell populations with known differences
  • predetermined criteria such as cancer cells and cancer cells, or other cell populations with known differences
  • the difference in SNP types and distribution among different taxa can effectively identify genes closely related to cell mutations, and further analyze the function of genes to determine the functions closely related to cell mutations, thereby determining cell mutations.
  • a specific state of an organism such as a human such as a marker of a disease, including a genetic marker and a functional marker.
  • the method which can be used for cell sorting is not particularly limited, and may be based on a clinical classification or a single cell sorting method as described above. It should be noted that the term "sub-population" is used herein to distinguish it from the "group" in the single-cell classification method. Without affecting the understanding, in this paper, sometimes the sub-population will be directly "called" group,,.
  • the present invention provides a gene screening device.
  • the device The apparatus includes: a calculation unit adapted to obtain a classified subpopulation according to a cell classification result, and calculate a statistic of all SNP sites of each gene in each type of population, optionally making a difference in the obtained statistic
  • the test obtains the test value; the sorting unit sorts all the genes by statistic or test value, and selects the gene with the highest statistic or test value as the gene related to the cell mutation.
  • the gene screening device of the present invention further comprises a cell sorting device which is the single cell sorting device described above to classify cells to obtain a classified population.
  • the present invention proposes a single cell classification method and screening method, and an apparatus for carrying out the method.
  • the reads (reading) results obtained by sequencing each single cell sample are compared with the reference genome sequence, and the alignment results are subjected to data filtering;
  • the uniform genotype of each single cell sample was determined based on the filtered data, and the consistent genotypes of all single cell samples were saved as SNP data sets;
  • a genotype file of the locus corresponding to the position of the reference genome SNP data set is extracted from the saved SNP data set; the cell mutation SNP site is selected, and the cells are classified according to the genotype file of the cell mutation SNP site.
  • a single cell sorting apparatus includes:
  • the data filtering module compares the reads obtained by sequencing each single cell sample with the reference genome sequence, and performs data filtering on the comparison result;
  • a genotyping module that determines a consistent genotype for each single cell sample based on the filtered data and saves the consistent genotype of all single cell samples as a SNP data set;
  • a genotype file extraction module that extracts a genotype file of a locus corresponding to a position of a reference genome SNP data set from a saved SNP data set;
  • the classification module selects the cell mutation SNP site and classifies the cells according to the genotype file of the cell population mutant SNP.
  • the classified population is obtained according to the cell classification result, the statistics of all SNP sites of each gene in each group are calculated, and the statistics are accumulated;
  • the predicted genes are sorted by statistics or test values, and the genes with the highest statistic or test value are selected.
  • the calculation unit obtains the classified population according to the cell classification result, calculates the statistics of all SNP sites of each gene in each group, and accumulates the statistics; performs a difference test on the obtained statistics to obtain the test value;
  • the sorting unit is coupled to the obtaining unit and the calculating unit, and sorts the predicted genes by statistic or test value, and selects the gene with the highest statistic or test value.
  • the present invention adopts next-generation sequencing technology (NGS), analyzes and studies single-cell genomes through bioinformatics methods, and collects sub-populations of cells (or 4) for in-depth follow-up research.
  • NGS next-generation sequencing technology
  • avoiding the labeling operation of cells effectively solves the problem that some cell subpopulations in the traditional single cell classification method do not have corresponding specific markers and cannot be classified; on the other hand, comprehensively and completely analyze single cells.
  • the genetic variation information of the genome and the setting of multiple control samples greatly improve the accuracy of cell subgroup classification.
  • the present invention also proposes a single cell screening method capable of screening out a subpopulation of cells (or children) to improve the accuracy of cell screening.
  • FIG. 2 is a schematic diagram of a unique contiguous reference genome of a prior art
  • FIG. 3 is a flow chart of a method for single cell classification and screening according to the present invention
  • Figure 4 is a N-J relationship tree of the renal cancer exome sequencing of the present invention.
  • Figure 5 is a maximum likelihood relationship tree of the renal cancer exome sequencing of the present invention.
  • Figure 6 is a diagram showing the result of PCA sequencing of the renal cancer exome of the present invention, wherein the abscissa represents the first principal component vector, and the ordinate represents the second principal component vector;
  • Figure 7 is a diagram showing the result of PCA sequencing of the renal cancer exome of the present invention, wherein the abscissa represents the first principal component vector, and the ordinate represents the third principal component vector;
  • Figure 8 is a diagram showing the result of PCA sequencing of the renal cancer exome of the present invention, wherein the abscissa represents the first principal component vector, and the ordinate represents the fourth principal component vector;
  • FIG. 9 is a diagram showing the results of the sequence of the renal cancer exome sequencing of the present invention, wherein "series 1" represents a cancer cell population, “series 2,” represents a paracancerous cell population;
  • Figure 10 is a schematic diagram showing the classification relationship between 53 cancer cells and 8 normal cells of the present invention.
  • FIG. 11 is a schematic diagram of clustering of cancer cells and normal cells of the present invention, wherein the abscissa represents a first principal component vector, and the ordinate represents a second principal component vector;
  • Figure 12 is a schematic view of a single cell sorting device of the present invention.
  • FIG. 13 is a schematic diagram of a screening module in the single cell sorting apparatus of the present invention. Detailed description of the invention
  • the present invention adopts next-generation sequencing technology (NGS), analyzes and studies single-cell genomes through bioinformatics methods, and screens and collects sub-populations of cells (or children) for further follow-up research.
  • NGS next-generation sequencing technology
  • clinical diagnosis and treatment such as prenatal diagnosis, preimplantation genetic diagnosis, individualized treatment, multi-point mapping, sperm and egg typing, genetic diagnosis, tumors (such as lymphoma) , leukemia), medical research (such as autism, neurological diseases and autoimmune diseases, genomic variation studies, stem cell research, drug resistance studies, drug target gene screening, transcriptome detection) , cell model studies, population identification, etc.), archaeological studies, forensic testing.
  • the single cell sample referred to in the present invention includes a nucleic acid (genomic DNA or RNA such as non-coding RNA, mRNA); the single cell is derived from an organism and is prepared by a conventional method.
  • DNA or RNA can be extracted or amplified from bacteria, protozoa, fungi, viruses, and single cells of higher organisms/high animals and plants, such as mammals, particularly humans.
  • Single cells can be cultured in vitro, or from clinical samples (including plasma, serum, spinal fluid, bone marrow, lymph, ascites, pleural effusion, oral fluid, skin tissue, respiratory tract, digestive tract, reproductive tract, urinary tract, tears, Directly isolated from saliva, blood cells, stem cells, tumors, fetal cells can be derived from embryos (such as one or several embryos/embryos) or maternal blood, and can come from living or dead organisms. Samples included single cell suspension, paraffin-embedded tissue sections, and biopsy tissue.
  • Samples can reflect cell-specific states such as cell proliferation, cell differentiation, apoptosis/death, disease status, external stimuli, and developmental stages.
  • single cell samples are derived from mammals, including preimplantation embryos, stem cells, suspected cancer cells, pathogenic organisms, and even from crime scenes.
  • analysis of human blastomeres can detect genetic defects in the fetal genome.
  • Figure 3 shows the flow from (7).
  • Lysis of cells The individual cells isolated are lysed by a detergent method, a boiling method, an alkali denaturation method, a lysozyme method, an organic solvent method, etc., to obtain intact cell genomic DNA.
  • PCR-based amplification such as DOP-PCR, PEP-PCR, T-PC
  • linear DNA amplification such as OmniPlex WGA, multiple displacement amplification (MDA).
  • Single-cell whole-genome amplification is performed to achieve the DN A starting amount required by next-generation sequencing technology.
  • Exo me sequence capture technology for DNA library construction. After quality inspection, single-cell genome sequencing is performed using next-generation sequencing technology, such as Illumina HiSeq 2000 sequencing system, Illumina Genome Analyzer II sequencing system, AB SOLiDTM 4.0 sequencing system, Roche GS FLX Titanium System, etc.
  • next-generation sequencing technology such as Illumina HiSeq 2000 sequencing system, Illumina Genome Analyzer II sequencing system, AB SOLiDTM 4.0 sequencing system, Roche GS FLX Titanium System, etc.
  • Bowtie is aligned with reference genome sequences (eg, human genomes HG18, HG19) and optimal parameters are set based on specific data to accurately locate Reads.
  • reference genome sequences eg, human genomes HG18, HG19
  • Sequencing depth refers to the average depth at which a genome is measured, calculated by dividing the total number of bases sequenced by the size of the genome.
  • Sequencing coverage refers to the approximate proportion of the portion of a genome that is measured. The calculation is the site to which the genome is covered divided by the effective length of the genome.
  • Sequencing depth and coverage were used to assess whether the amount of data used for bioinformatics analysis was sufficient and whether sequencing was uniform.
  • a pair of reads are randomly selected when the sequences of the pairs of reads are completely identical, and other repeated double-end reads in the plurality of pairs of reads are removed; and/or the unpaired comparison reference genome is removed.
  • the repeated p end-end reads in each DNA library are selected, for example, repeated double-end reads due to excessive PCR amplification times, of course, not limited to PCR amplification, but also other expansions.
  • the addition mode results in multiple pairs of repeated double end reads.
  • the three pairs of reads A, B, and C sequences are identical and aligned exactly to the start and end positions on the genome, and the start and end positions are identical and have duplicate double ends. In this case, only one pair of reads is randomly reserved, Other duplicate reads are removed.
  • exome sequencing of the human genome is not limited to this, for example, sequencing of other mammals, and the sequencing method is not limited to exome sequencing, such as whole genome sequencing.
  • the reads from the exome sequencing should be mostly aligned to the human reference genome.
  • only the reads of the unique reference genome are selected (ie reads with a hit number of 1), which greatly reduces the impact of the error comparison.
  • Reads D aligns multiple positions on the reference genome, while Reads E only compares the position to the top. Since the exome is not a repeating region in the genome, Reads D is removed directly.
  • genotyping software such as SOAPsnp, SAMtools Pileup pileup
  • the target region is the region in which the exon is located. Specific areas that require sequencing and bioinformatics analysis are typically indicated:
  • the present invention synthesizes the probability files for all cell genotypes, integrates all cell data using the maximum likelihood approach, and produces pseudo-containing all sites for each site.
  • Genome Pseudo-genome. Select the genotype with the highest probability as the consistent genotype of each cell, and detect high-quality SNPs by information such as genotype and sequencing depth. After obtaining the consistent sequence of samples, the results are saved as SNP datasets in the population SNPs format.
  • the genotype of the corresponding site of each cell is extracted from the genotype possibility file, and the genotype file of the population SNP at the corresponding site of each cell is obtained.
  • the locus refers to the location ( osition ) where the SNP is located.
  • the present invention is mainly aimed at finding a difference site between cells, and therefore it is necessary to select a site related to cell mutation. First, the SNP sites between non-populations are removed.
  • SNP loci between non-populations All individual base types are the same, and the relative reference sequence is a SNP, then the site is a non-community SNP locus.
  • the reference sequence is A, and all individuals are heterozygous base type R at this position, which is a non-community SNP site.
  • the published SNP locus is removed.
  • the normal human SNP locus is removed, that is, the dbSNP of the human genome HG 18, the SNP of Yanhuang No. 1 and the SNP of the thousand human genome are removed.
  • the obtained SNP site is a cell mutation SNP site.
  • the methods commonly used in the analysis of bioinformatics are used, such as the adjacent NJ method, MEGA software, principal component analysis (Princi al Com onents Analysis, PCA X group structure (Group Structure) Sorting the cells, etc. At least one of the above methods may be employed in performing cell sorting. As an embodiment of the present invention, all of the above methods are employed, and when the classification results of the various methods are identical, the final result is confirmed. Cell classification results.
  • D represents the difference distance between cells and cells, where "is the total number of SNPs , and ⁇ indicates the degree of difference between cells and cells at the site, where
  • genotypes are identical, such as position, cell: , ':
  • Genotypes are completely different, such as position, cells, ': C
  • the genotypes are different, such as position, cells: ⁇ , j': M Since the human genome is diploid, A indicates that both loci in the allele are A, M is a heterozygous locus, ie, A and C combination.
  • the Fneighbor program http: ⁇ emboss.bioinformatics.nl/cgi-bin/ emboss help/fneighbor ), which obtains its evolutionary tree by adjacency (NJ).
  • the genotype files of all SNP loci in each cell are composed into a sequence, one cell corresponds to one sequence, and as an input file of MEGA, MEGA is based on each cell sequence.
  • the difference, in which the software has three methods (M aximum likelihood, Least Squares and Maximum p ars imony ) constitute a relationship tree.
  • Principal Composits Analysis is a technique for binarizing data sets and is a linear transformation. This transformation transforms the data into a new coordinate system such that the first large variable of any data projection is on the first coordinate (called the first principal component) and the second largest variable is in the second coordinate (second On the main component), and so on.
  • Principal component analysis is often used to reduce the dimensionality of a data set while preserving the feature variables that contribute the most to the data set. This is achieved by preserving the low-order principal components and ignoring the higher-order principal components. This is because low-order components tend to preserve the most important aspects of the data set.
  • the SNP data is first converted into a digital matrix (consistent with the reference sequence is 0, the opposite is 2, the base is 1) and make uniformity. Then construct a linear vector equation by the method described above.
  • the Structure software ( htt ://prit ch.bsd.uchicago.edu/software/structure2 1. html ) is based on the genotyping data of the SNP locus to infer whether there are different groups and determine the group to which each individual belongs. According to the software description, the genotype file conversion format of the population SNP is used as the Structure input file and up to 50,000 simulations are used in the hybrid model. The probability of each cell belonging to each group is calculated assuming that multiple groups exist.
  • single cells can be further screened as follows:
  • the classification of single cells is realized, and the information of each cell sample is extracted, and the cells with objections are excluded, such as unclear classification or obvious outlier samples.
  • the genomes of these groups can be compared through a series of statistics and tests to find regions or genes with significant differences, and genes with higher correlation coefficients can be screened.
  • the annotation file corresponding to HG18 was downloaded from the Human Genome Database to obtain the starting and ending positions of more than 30,000 genes in the human genome that have been predicted.
  • Each gene referred to herein refers to a gene in a predicted genome.
  • the formula for calculating statistics is mainly used to measure the polymorphism of a group.
  • a ' b refers to the number of samples of two bases in a group.
  • the formula can be:
  • C 2 h can also perform a difference test on the obtained statistics to obtain the test value.
  • the test value used is at least one of the following: Lod, Fst, Pbs.
  • the above three test values may be used, and when the above three test values are identical, the result is the final test value.
  • the last method may be employed, and when the sort result obtained based on the statistic coincides with the sort result obtained based on the test value, it is used as the finally selected gene.
  • FIG. 12 it is a schematic diagram of a single cell sorting device of the present invention.
  • the device includes:
  • the data filtering module compares the read results obtained by sequencing each single cell sample with the reference genome sequence, and performs data filtering on the comparison result;
  • a genotyping module coupled with the data filtering module, determines a genotype of each single cell sample based on the filtered data, and saves the consistent genotype of all single cell samples as a SNP data set;
  • a genotype file extraction module coupled with the genotype determination module to extract a genotype file of a locus corresponding to a position of the reference genomic SNP data set from the saved SNP data set;
  • the classification module is coupled with the genotype file extraction module to select a cell mutation SNP site, and classifies the cells according to the genotype file of the cell population mutation SNP, and the classification method adopted includes at least one of the following: contiguous NJ method framing, MEGA Software, PCA principal component analysis, Structure group structure.
  • the single cell classification device further includes:
  • a screening module coupled with the classification module, for obtaining a start and stop position of the gene in the predicted genome; obtaining the classified population according to the classification result, calculating statistics of all SNP sites of each gene in each group, and accumulating statistics; The obtained statistic is used for the difference test to obtain the test value; the predicted genes are sorted by statistic or test value, and the gene with the highest statistic or test value is selected.
  • the screening module may further include the following units, as shown in FIG. 13, wherein:
  • the calculation unit obtains the classified population according to the cell classification result, calculates the statistics of all SNP sites of each gene in each group, and accumulates the statistics; performs a difference test on the obtained statistics to obtain the test value;
  • a sorting unit coupled to the obtaining unit and the calculating unit, sorting the predicted genes by statistics or test values, Screen out the genes with the highest statistic or test value.
  • RC-1 to RC-15 represent renal cell single cells, a total of 15 single cell samples; RN-1 to RN-5 represent adjacent cancer cells; RN-T indicates normal tissue direct extraction DNA sequencing, used as a control Data analysis and evaluation. Main Single cancer cells were used as control samples. There are also cases in which the adjacent cancer cells and normal tissues are simultaneously used as control samples, such as the above two control samples when the site of loss of heterozygosity is removed.
  • the sequencing depth was increased in the exome sequencing.
  • the average depth of the exon target region (Mean Depth) was around 30x, the target region coverage reached 80-96%.
  • a site with four reads support can determine the accuracy of the base of the site to 99%, and the statistics show that the site with a sequencing depth of 4 accounts for 60-90 of the exon region. % , indicating that the data from the exome sequencing is better than the data from the whole genome sequencing.
  • Table 2 Coverage and depth data of exon sequencing of renal cell single-cell samples
  • the depth of whole genome sequencing is low and the subsequent analysis cannot be performed, and the depth of sequencing by the exome is high.
  • the following is mainly based on the data obtained by sequencing of exomes.
  • the repeated double-end reads caused by too many PCR amplifications in each DNA library are selected.
  • a pair of reads are randomly selected to remove other reads.
  • the three pairs of reads A, B, and C sequences in Figure 1 are identical and aligned to the start and end positions on the genome. In this case, only one pair of reads is randomly reserved.
  • the human exon region is not likely to have multiple copies on the genome, that is, it is impossible to repeat the sequence, so the reads from the exome sequencing should be mostly unique to the human reference genome. .
  • the reads of the unique reference genome are selected (ie reads with a hit number of 1), which greatly reduces the impact of the error comparison.
  • Reads D aligns multiple locations on the reference genome, while Reads E only compares the position—only because the exon is not a repeat region in the genome, so Reads D is removed directly. 1-4, determination of individual genotype
  • the filtered data in steps 1-3 is input into the SOAPsnp software, and the filtered data is compared with the reference genome by the SOAPsrp software, and the comparison result is a CNS file.
  • the present invention synthesizes the possibility of all cell genotypes, integrates all cell data using the maximum likelihood method, and produces a pseudogenome containing each site for all samples (Pseudo- Genome ) facedSelect the most probable genotype as the consistent genotype of each cell, and detect high quality SNPs by genotype and depth information. After obtaining the consistent sequence of the sample, the results are saved as a population SNPs format.
  • the corresponding sites of each cell genotype were extracted from the CNS file, and the genotype files of the population SNPs at the corresponding sites of each cell were obtained.
  • the file format is shown in Table 3.
  • SNP position indicates the position of the SNP site on the chromosome
  • base type corresponds to the base type of a certain cell genome at this position
  • site with a depth of 0 is represented by "-" (ie, data loss) Site).
  • sample ID corresponds to 21 single cell or tissue DNA samples.
  • Table 3 The genotype file format of the SNP in the corresponding sites of each cell
  • the number of SNPs in the 21 single-cell or tissue DNA samples relative to the human genome HG18 in the target region was 93,957.
  • the combination of heterozygous loci is represented by the following letters:
  • the present invention is mainly aimed at finding a difference site between cells, and therefore it is necessary to select a site related to cell mutation.
  • the ninth single-cell sample RC-9 is likely to be originally a heterozygous site at this site, but since it is only amplified to one base type, it is judged to be a homozygous site.
  • the normal tissue directly extracted DNA-sequenced RN-T (ie, the last column) must be a heterozygous site. Because this sample is not amplified by WGA, no loss of heterozygosity is possible.
  • the number of heterozygous sites and data loss samples in the other 20 single-cell samples is greater than or equal to 18. That is, only one single cell sample is allowed to be heterozygous at the same site. Because three or more samples are simultaneously at the same site, the probability of heterozygous loss is extremely low.
  • the number of SNPs is 89,478.
  • the SNP site of the normal human being released must be removed, that is, the dbSNP of the human genome HG18, the SNP of Yanhuang No. 1 and the SNP of the thousand human genome are removed, and There were 50,524 SNP loci associated with cell mutation.
  • the cells were classified by the methods commonly used in bioinformatics population analysis. Classification is determined by the branching and clustering of the phylogenetic tree. As shown in Figure 4, RC and RN are clearly clustered into two separate parts in the phylogenetic tree, so they are divided into two categories.
  • cells can be classified according to the phylogenetic tree.
  • Figure 5 is a relational tree constructed by the maximum likelihood method, which classifies cells according to the relationship tree.
  • the abscissa indicates the sample name
  • the ordinate indicates the probability of the group to which each sample belongs, and the cells are classified according to the probability.
  • 20 single cells can be roughly divided into two groups. Renal cancer exome sequencing The results of the Structure are shown in Figure 9.
  • each cell sample information was extracted, and the dissenting cells (unclear or apparent outlier samples) were excluded. From the results of the above various group structure analysis, the sampling is normal and the classification is reasonable. These 20 single-cell samples can be roughly divided into two groups, namely cancer cell population (15 RC) and paracancerous cell population (5 RN), of which RC-1, RC-6, and RC-7 are cancer cells. A subgroup.
  • the cell sample information refers to which of the single cells to be analyzed are cancer cells and which are cancerous cells (determined by sampling).
  • the cell sample information is only for reference and needs to be analyzed together with the clustering results. If the cell sample information at the time of sampling is considered to be a cancer cell and a cancerous cell, and clustering is divided into two populations, the results correspond to each other; if the cell sample information at the time of sampling is inconsistent with the clustering result, the clustering result is quasi.
  • 1-11.1 Download the annotation file corresponding to HG18 from the Human Genome Database to obtain the starting and ending positions of more than 30,000 genes in the currently predicted human genome.
  • the formula for calculating statistics is mainly used to measure the polymorphism of a group.
  • a and b refer to the number of samples of two bases in a group.
  • the formula can be:
  • test value used is at least one of the following:
  • the test values are ⁇ , Fst, Pbs, and the above three test values are used in this embodiment.
  • the individual test values and the calculation process are specifically described below.
  • ⁇ is as follows: If ⁇ ⁇ and the difference is small, that is, the gene does not differ much in the two populations, Lod (which can deviate significantly from 0) can be preliminarily thought to be an important gene that causes differentiation between the two populations. , respectively, the value of more than 30,000 genes in HG18 is £c3 ⁇ 4 g ⁇ , and then sorted from large to small to select the top ranked genes.
  • FsT Fitification index
  • ⁇ ⁇ ⁇ means that two samples are taken from a group (Within) to form a pair, and the difference of the SNP genotypes of the pair of samples is calculated, so that the difference of SNP genotypes of all pairs of samples can be calculated, and finally the average is Nwi in . If there are two groups, the two groups calculate 11 and then accumulate.
  • X is the frequency of the sub-allocation base (second base) of the SNP locus i in the population j; but the physical position of the SNP locus on the chromosome in the population; then the population is used for analysis.
  • the variables take RC and RN, and the variables are substituted with the last determined SNP position.
  • the value of each gene is calculated in units of genes, and then the values of more than 30,000 genes in HG18 are sorted to select the top ranked genes.
  • T RC - B - ⁇ g( ⁇ - Fst RC _ B )
  • the embodiment obtains important genes according to ⁇ and corpse & screening, and performs functional analysis separately. It is judged whether these genes are affected in certain pathways and thus related to the mechanism of renal cancer.
  • Example 2 Classification and screening of leukemia single cells
  • a 30x depth exome sequence was sequenced for each cancer single cell, and the resulting reads were aligned with the reference genomic sequence (Human Genome HG18) using the SOAPaligner 2.0 alignment software. Since the human SNP is in the order of two thousandths and the reading length of Reads is about lOObp, we set up at most 2 error comparisons (mismacth) for each Reads in SOAP comparison, and Gap is not allowed to ensure the comparison to the reference. The accuracy of Reads on the genome.
  • the consensus genotype of each locus in the Soap snp software has a mass value of not less than 20, and the rank test p value is not less than 1%; for heterozygous variant SNPs: the locus genotype differs from the reference genome, preferably the base type
  • the sequencing quality value is not less than 20, the sequencing depth is not less than 6, the sub-good base type sequencing quality value is not less than 20, the sequencing depth is not less than 2, and the ratio of the two genotype sequencing depths is in the range of 0.2-5.
  • Fig. 10 is a schematic diagram showing the classification relationship between 53 cancer cells and 8 normal cells of the present invention, wherein ET-T 1 represents cancer tissue and NC-T 1 represents normal tissue.
  • FIG. 11 is a schematic diagram showing clustering of cancer cells and normal cells of the present invention, wherein LC represents cancer cells and LN represents normal cells.
  • the cell sample information is extracted, and the dissenting cells (unclear classification, or outlier samples) are excluded. From the above group structure description, the product is normal and the classification is reasonable.
  • Normal cells have a consistent homozygous genotype, and there are two or more heterozygous or homozygous mutations in cancer cells, and there is no third homozygous genotype and heterozygosity that is inconsistent with the two homozygous genotypes. genotype. For example, if the normal cell genotype is A, or the mutation type is A->C, then only three genotypes, namely, C and M, can be found in cancer cells, and the number of C and M is not less than 2. This site is called High confidence somatic mutation (HCSM).
  • HCSM High confidence somatic mutation
  • the position information of the gene mutation site and the number of non-synonymous mutation sites in each gene are used as the standard for gene function enrichment, and the Webgestalt online analysis tool is used.
  • the mutation affects the function and pathway of the gene (http://bioinfo.vanderbilt.edu/webgestalt/option.php), and it was found that the mutation was mainly concentrated in the following eight types of functional genes.
  • ECM -receptor interaction ECM receptor interaction
  • Viral myocarditis (viral '
  • MAPK signaling athwi ⁇ MAPK signaling pathway
  • Pantothenate and CoA biosynthesis biosynthesis of pantothenic acid and coenzyme A
  • CAMs Cell adhesion molecules

Landscapes

  • Life Sciences & Earth Sciences (AREA)
  • Chemical & Material Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Health & Medical Sciences (AREA)
  • Proteomics, Peptides & Aminoacids (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Organic Chemistry (AREA)
  • Analytical Chemistry (AREA)
  • General Health & Medical Sciences (AREA)
  • Biotechnology (AREA)
  • Biophysics (AREA)
  • Evolutionary Biology (AREA)
  • Wood Science & Technology (AREA)
  • Zoology (AREA)
  • Molecular Biology (AREA)
  • Theoretical Computer Science (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Medical Informatics (AREA)
  • Genetics & Genomics (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Biochemistry (AREA)
  • Immunology (AREA)
  • General Engineering & Computer Science (AREA)
  • Microbiology (AREA)
  • Measuring Or Testing Involving Enzymes Or Micro-Organisms (AREA)

Abstract

提供了单细胞分类方法、基因筛选方法以及实施该方法的装置。其中，单细胞分类方法包括以下步骤：对来自相同群体的多个单细胞样本的全基因组分别进行测序，以便获得来自每个单细胞样本的reads；分别将来自每个单细胞样本的reads与参考基因组序列进行比对，并对所述reads进行数据过滤；基于经过过滤的reads，确定每个单细胞样本的一致基因型，其中，所有单细胞样本的一致基因型构成所述群体的SNP数据集；针对所述每个单细胞，基于所述群体的SNP数据集，确定在与参考基因组SNP数据集的位置相应的位点，每个细胞所对应的基因型；以及选择与细胞突变相关的SNP位点，并且基于所述单细胞在该位点的基因型，对所述单细胞进行分类。

Description

单细胞分类方法、 基因筛选方法及其装置 优先权信息
本申请请求 201 1 年 8 月 25 日向中国国家知识产权局提交的、 专利申请号为 2011 10245356.8的专利申请的优先权和权益， 并且通过参照将其全文并入此处。 技术领域
本发明涉及生物信息学，尤其涉及单细胞分类和基因筛选方法以及用于所述方法的 装置。 背景技术
不同个体之间， 个体的不同组织之间， 甚至同一组织的不同部位在基因表达、 拷贝 数变异、表观遗传等方面都存在显著差异。 细胞之间也存在异质性， 即使是体外培养遗 传背景完全相同的细胞群体。对于干细胞或前体细胞，因为任何状态改变都是可遗传的， 细胞异质性尤为明显。 为了更好地研究细胞生物学， 揭示细胞异质性的规律， 非常需要 开发应用于单个细胞研究的技术方法， 因此有学者提出"单细胞分析（ SCA ) "概念， 从 "组学 （Omics ) "角度进行阐述。 单细胞分类和筛选为单细胞分析提供了重要基础。
单细胞分类可以有效应用于各种干细胞分化过程的研究中， 如肿瘤干细胞、胚胎干 细胞的定向分化、 造血干细胞的研究中， 需要筛选不同分化阶段的干细胞， 进行各种干 细胞的检测。在耐药性研究中， 需要对给药不同时期的细胞进行精确分类， 从而进一步 分析该细胞亚群的耐药性和耐药基因，例如可进行癌症病人的多药耐药性及多药耐药基 因与药物滥用、 药物耐受、 药物依赖的关系的研究。 同样地， 在药物靶点基因的筛选中， 由于药物与细胞， 特别是敏感细胞相互作用，将引起细胞外部形态及内部正常代谢过程 的一系列变化， 因此筛选出敏感细胞是关键的第一步， 为后期精确定位药物靶点基因提 供重要基础。 单细胞分类和筛选应用于建立药效筛选模型， 为药物设计、靶点的选择和 用药方案的确定提供理论依据， 同时使药物筛选有了更高的特异性。
目前， 常用的筛选单细胞方法多为物理机械、 化学或生物的方法， 如流式细胞仪、 磁性细胞分选仪等方法。一方面， 这些技术采用表面活性剂、 荧光染料、 抗原抗体， 细 胞毒性大， 只能对特异标记的或非特异标记的单细胞悬液进行分选，前期样本制备过程 繁瑣， 且目前对众多荧光探针、 单抗（包括细胞表面 CD分子） 的特异性争论较多， 许 多细胞亚群并无对应的特异性标记物 /特异性抗原； 另一方面， 这些技术采用生物学、 免疫学、 化学方法， 通过表型测定（包括细胞大小、 细胞粒度、 细胞表面积、 核浆比例 等） ， 进行统计学分析， 对于亚群分类、 筛选和检测的灵敏度低， 缺乏有效的准确性评 估。 发明内容
本发明旨在至少解决现有技术中存在的技术问题之一。
在本发明中， 除另有说明， 否则本文中使用的科学和技术术语具有本领域技术人员所 通常理解的含义。 同时， 为了更好地理解本发明， 下面提供相关术语的定义和解释。
术语 "基因型的可能性文件" , 是指利用 SNP检测软件，设置先验概率参数利用贝叶斯 公式计算出的样本目标区域可能的基因型的后验概率的数值集合； 当利用的 SNP检测软件 是 SOAPsrp时， 生成的 "基因型的可能性文件" 即为 CNS文件。
如本文使用的， "基因型文件" 是指选择上述 "基因型的可能性文件" 中概率最大的基 因型作为每个细胞的一致基因型后， 根据参考基因组 SNP数据集位置信息， 提取每个细胞 基因型的相应位点， 获得的群体 SNP在各细胞相应位点的基因型集合。
在本发明的一个方面， 本发明提出了一种单细胞分类方法。 根据本发明的实施例， 本 发明的单细胞分类方法包括： 对来自相同群体的多个单细胞样本的全基因组分别进行测序， 以便获得来自每个单细胞样本的 reads; 分别将来自每个单细胞样本的 reads与参考基因组 序列进行比对， 并对所述 reads进行数据过滤； 基于经过过滤的 reads , 确定每个单细胞样 本的一致基因型， 其中， 所有单细胞样本的一致基因型构成所述群体的 SNP数据集； 针对 所述每个单细胞， 基于所述群体的 SNP数据集， 确定在与参考基因组 SNP数据集的位置相 应的位点， 每个细胞所对应的基因型； 以及选择与细胞突变相关的 SNP位点， 并且基于所 述单细胞在该位点的基因型， 对所述单细胞进行分类。 由此， 根据本发明的实施例， 可以 采用新一代测序技术（NGS ), 通过生物信息学方法， 对单细胞基因组进行分析和研究， 收 集细胞亚群（或微粒） 进行深入的后续研究。 一方面， 避免对细胞进行标记的操作， 有效 解决了传统的单细胞分类方法中某些细胞亚群无对应的特异性标记物而无法进行分类的问 题； 另一方面， 全面完整地分析单细胞基因组的遗传变异信息， 并可以设置多个对照样本， 大大提高了细胞亚群分类的准确性。
根据本发明的实施例， 上述单细胞分类方法还可以具有下列附加技术特征：
在本发明的一个实施例中， 利用第二代或者第三代测序平台进行所述测序， 所述数据 过滤的标准为： 当存在多对重复双末端 reads时， 在多对 reads的序列完全一致时随机选择 一对 reads, 并将所述多对 reads中其他重复双末端 reads去除； 和 /或去除未唯一比对上所述 参考基因组序列的 reads。
在本发明的一个实施例中， 基于经过过滤的 reads , 确定每个单细胞的一致基因型进一 步包括： 基于所述经过过滤的 reads,确定每个单细胞样本在目标区域基因型的可能性； 基 于所有单细胞样本在目标区域基因型的可能性， 确定包含所有样本每个位点的伪基因组； 以及从所述伪基因组中选择概率最大的基因型作为每个单细胞样本的一致基因型。
在本发明的一个实施例中，选择与细胞突变相关的 SNP位点进一步从所述群体 SNP数 据集中去除下列的至少之一： 非群体间的 SNP位点、 存在杂合丢失的位点、 以及已经公布 的 SNP位点。
在本发明的一个实施例中， 所述多个单细胞样本至少之一的全基因组在进行测序之前 被进行全基因组扩增处理， 其中， 去除存在杂合丢失的位点进一步包括去除符合下列条件 的位点： 针对未经全基因组扩增的样本， 测序结果为杂合位点； 并且在经过全基因组扩增 的样本中， 在同一位点， 杂合位点和数据丢失的样本数目大于等于经过全基因组扩增的样 本的数目减去 3。
在本发明的一个实施例中， 针对所述每个单细胞， 基于所述群体的 SNP数据集， 确定 在与参考基因组 SNP数据集的位置相应的位点， 每个细胞所对应的基因型进一步包括按照 下列标准对所述 SNP数据集进行筛选： 每个位点的一致性基因型的质量值不小于 20, 秩检 验 p值不小于 1%; 以及对于杂合变异的 SNPs: 最好碱基型测序质量值不小于 20、 测序深 度不小于 6, 次好碱基型测序质量值不小于 20, 测序深度不小于 2, 且两基因型测序深度的 比值在 0.2-5范围内。
在本发明的一个实施例中， 对细胞进行分类之后， 还包括以下步骤： 提取每个细胞样 本信息， 剔除存在异议的细胞。
在本发明的一个实施例中， 对所述单细胞进行分类之后， 进一步包括： 基于分类结果 确定已分类群体， 并计算每类群体中每个基因所有 SNP位点的统计量， 任选地对所获得的 统计量作差异检验， 获得检验值； 选择统计量或检验值最高的基因作为与细胞突变相关的 基因。
在本发明的又一方面， 本发明提出了一种单细胞分类装置。 根据本发明的实施例， 该 单细胞分类装置包括： 数据过滤模块， 所述数据过滤模块适于将分别将来自每个单细胞样 本的 reads与参考基因组序列进行比对， 并对所述 reads进行数据过滤， 其中， 所述每个单 细胞样本的 reads是通过对多个单细胞样本的全基因组分别进行测序获得的； 基因型确定模 块， 所述基因型确定模块适于基于经过过滤的 reads, 确定每个单细胞样本的一致基因型， 其中， 所有单细胞样本的一致基因型构成所述群体的 SNP数据集； 基因型文件提取模块， 所述基因型文件提取模块适于针对所述每个单细胞， 基于所述群体的 SNP数据集， 确定在 与参考基因组 SNP数据集的位置相应的位点， 每个细胞所对应的基因型； 分类模块， 所述 分类模块适于基于预先选择的与细胞突变相关的 SNP位点， 并且基于所述单细胞在该位点 的基因型， 对所述单细胞进行分类。 利用该装置能够有效地实施前述单细胞分类方法。 由 此， 根据本发明的实施例， 可以采用新一代测序技术（NGS ), 通过生物信息学方法， 对单 细胞基因组进行分析和研究， 收集细胞亚群（或微粒） 进行深入的后续研究。 一方面， 避 免对细胞进行标记的操作， 有效解决了传统的单细胞分类方法中某些细胞亚群无对应的特 异性标记物而无法进行分类的问题； 另一方面， 全面完整地分析单细胞基因组的遗传变异 信息， 并可以设置多个对照样本， 大大提高了细胞亚群分类的准确性。
根据本发明的实施例， 该单细胞分类装置还可以具有下列附加技术特征：
在本发明的一个实施例中， 所述数据过滤模块适于基于下列标准进行数据过滤： 当存 在多对重复双末端 reads时， 在多对 reads的序列完全一致时随机选择一对 reads, 并将所述 多对 reads 中其他重复双末端 reads 去除； 和 /或去除未唯一比对上所述参考基因组序列的 reads 在本发明的一个实施例中， 所述基因型确定模块适于通过下列确定所述每个单细胞的 一致基因型： 基于所述经过过滤的 reads,确定每个单细胞样本在目标区域基因型的可能性； 基于所有单细胞样本在目标区域基因型的可能性， 确定包含所有样本每个位点的伪基因组； 以及从所述伪基因组中选择概率最大的基因型作为每个单细胞样本的一致基因型。
在本发明的一个实施例中， 分类模块适于从所述群体 SNP数据集中去除下列的至少之 一， 选择与细胞突变相关的 SNP位点： 非群体间的 SNP位点、 存在杂合丢失的位点、 以及 已经公布的 SNP位点。
在本发明的一个实施例中， 所述多个单细胞样本至少之一的全基因组在进行测序之前 被进行全基因组扩增处理， 其中， 所述分类模块适于去除符合下列条件的位点， 以便去除 存在杂合丢失的位点： 针对未经全基因组扩增的样本， 测序结果为杂合位点； 并且在经过 全基因组扩增的样本中， 在同一位点， 杂合位点和数据丢失的样本数目大于等于经过全基 因组扩增的样本的数目减去 3。
在本发明的一个实施例中， 所述基因型文件提取模块适于按照下列标准对所述 SNP数 据集进行筛选： 每个位点的一致性基因型的质量值不小于 20, 秩检验 p值不小于 1%; 以及 对于杂合变异的 SNPs: 最好碱基型测序质量值不小于 20、 测序深度不小于 6, 次好碱基型 测序质量值不小于 20, 测序深度不小于 2, 且两基因型测序深度的比值在 0.2-5范围内。
在本发明的一个实施例中， 所述分类模块进一步适于提取每个细胞样本信息， 剔除存 在异议的细胞。
在本发明的一个实施例中， 进一步包括筛选模块， 所述筛选模块适于： 基于分类结果 确定已分类群体， 并计算每类群体中每个基因所有 SNP位点的统计量， 任选地对所获得的 统计量作差异检验， 获得检验值； 以及选择统计量或检验值最高的基因， 作为与细胞突变 相关的基因。
在本发明的再一方面， 本发明提出了一种基因筛选方法。 根据本发明的实施例， 该方 法包括以下步骤： 对细胞进行分类， 以便获得已分类亚群体， 并计算每类亚群体中每个基 因所有 SNP位点的统计量， 任选地对所获得的统计量作差异检验， 获得检验值； 选择统计 量或检验值最高的基因作为与细胞突变相关的基因。 通过预先对细胞进行分类， 例如可以 按照预定的标准， 诸如分为癌旁细胞和癌细胞， 或者其他具有已知区别的细胞群体， 通过 对每类群体中的 SNP位点进行统计分析，例如根据在不同分类群体中 SNP类型及分布的差 异， 可以有效地确定与细胞突变密切相关的基因， 并且进一步通过对基因的功能进行分析， 可以确定与细胞突变密切相关的功能， 由此可以确定细胞突变或者生物体例如人的特定状 态诸如疾病的标记物， 包括基因标记物和功能标记物。 才艮据本发明的实施例， 可以用于进 行细胞分类的方法并不受特别限制， 可以依据临床分类， 也可以为前面所描述的单细胞分 类方法。 需要说明的是， 在本文中所使用的术语 "亚群体" 是为了与单细胞分类方法中的 "群体"进行区分， 在不影响理解的前提下， 在本文中， 有时也直接将 "亚群体"称为 "群 体，，。
在本发明的再一方面， 本发明提出了一种基因筛选装置。 根据本发明的实施例， 该装 置包括： 计算单元， 所述计算单元适于根据细胞分类结果获得已分类亚群体， 并计算每类 群体中每个基因所有 SNP位点的统计量， 任选地对所获得的统计量作差异检验， 获得检验 值； 排序单元， 所述排序单元将所有基因按统计量或检验值进行排序， 筛选出统计量或检 验值最高的基因作为与细胞突变相关的基因。 利用该装置， 可以有效地实施前面所述的基 因筛选方法， 通过预先对细胞进行分类， 例如可以按照预定的标准， 诸如分为癌旁细胞和 癌细胞， 或者其他具有已知区别或者具有显著统计差异的细胞群体， 通过对每类群体中的 SNP位点进行统计分析， 例如根据在不同分类群体中 SNP的类型或分布的差异， 可以有效 地确定与细胞突变密切相关的基因， 并且进一步通过对基因的功能进行分析， 可以确定与 细胞突变密切相关的功能， 由此可以确定细胞突变或者生物体例如人的特定状态诸如疾病 的标记物， 包括基因标记物和功能标记物。 才艮据本发明的实施例， 细胞分类结果可以通过 前面所述的单细胞分类方法进行实施。 由此， 根据本发明的实施例， 本发明提出的基因筛 选装置进一步包括细胞分类装置， 该细胞分类装置为前面所述的单细胞分类装置， 以便对 细胞进行分类， 以便获得已分类群体。
由此， 根据本发明的实施例， 鉴于现有单细胞的分类和筛选方法存在的问题， 本发明 提出了单细胞分类方法和筛选方法， 以及实施所述方法的装置。
才艮据本发明的实施例的单细胞分类方法， 包括以下步骤：
将每个单细胞样本经测序得到的 reads (读段）结果与参考基因组序列进行比对， 并将 比对结果进行数据过滤；
根据过滤后的数据确定每个单细胞样本的一致基因型 （genotype ), 并将所有单细胞样 本的一致基因型保存为 SNP数据集；
从已保存的 SNP数据集提取与参考基因组 SNP数据集位置对应的位点的基因型文件； 挑选细胞突变 SNP位点， 并根据细胞突变 SNP位点的基因型文件， 对细胞进行分类。 根据本发明的实施例的单细胞分类装置， 包括：
数据过滤模块， 将每个单细胞样本经测序得到的 reads与参考基因组序列进行比对， 并 将比对结果进行数据过滤；
基因型确定模块， 根据过滤后的数据确定每个单细胞样本的一致基因型， 并将所有单 细胞样本的一致基因型保存为 SNP数据集；
基因型文件提取模块，从已保存的 SNP数据集提取与参考基因组 SNP数据集位置对应 的位点的基因型文件；
分类模块， 挑选细胞突变 SNP位点， 根据细胞群体突变 SNP的基因型文件， 对细胞进 行分类。
才艮据本发明的实施例的单细胞筛选方法， 包括以下步骤：
得到已预测基因组中基因的起止位置；
根据细胞分类结果获得已分类群体，计算每类群体中每个基因所有 SNP位点的统计量， 并累加统计量；
对获得的统计量作差异检验， 获得检验值； 将已预测基因按统计量或检验值进行排序， 筛选出统计量或检验值最高的基因。
根据本发明的实施例的单细胞筛选装置， 包括：
获取单元， 得到已预测基因组中基因的起止位置；
计算单元， 根据细胞分类结果获得已分类群体， 计算每类群体中每个基因所有 SNP位 点的统计量， 并累加统计量； 对获得的统计量作差异检验， 获得检验值；
排序单元， 耦合于获取单元和计算单元， 将已预测基因按统计量或检验值进行排序， 筛选出统计量或检验值最高的基因。
本发明采用新一代测序技术（NGS ), 通过生物信息学方法， 对单细胞基因组进行分析 和研究， 收集细胞亚群 （或 4 粒）进行深入的后续研究。 一方面， 避免对细胞进行标记的 操作 , 有效解决了传统的单细胞分类方法中某些细胞亚群无对应的特异性标记物而无法进 行分类的问题； 另一方面， 全面完整地分析单细胞基因组的遗传变异信息， 并设置多个对 照样本， 大大提高了细胞亚群分类的准确性。
本发明还提出单细胞筛选方法， 能够筛选出细胞亚群（或孩 ί粒）， 提高了细胞筛选的准 确性。
本发明的附加方面和优点将在下面的描述中部分给出，部分将从下面的描述中变得 明显， 或通过本发明的实践了解到。 附图说明
本发明的上述和 /或附加的方面和优点从结合下面附图对实施例的描述中将变得明 显和容易理解， 其中：
图 1为现有技术的重复片段 （Duplication Reads ) 示意图；
图 2为现有技术的唯一比对上参考基因组的片段（Unique mapped reads )的示意图； 图 3为本发明单细胞分类和筛选的方法流程图；
图 4为本发明肾癌外显子组测序的 N- J关系树；
图 5为本发明肾癌外显子组测序的最大似然关系树；
图 6为本发明肾癌外显子组测序 PCA结果图， 横坐标表示第一主成分向量， 纵坐 标表示第二主成分向量；
图 7为本发明肾癌外显子组测序 PCA结果图， 横坐标表示第一主成分向量， 纵坐 标表示第三主成分向量；
图 8为本发明肾癌外显子组测序 PCA结果图， 横坐标表示第一主成分向量， 纵坐 标表示第四主成分向量；
图 9为本发明肾癌外显子组测序 Structure结果图，其中 "系列 1 "表示癌细胞群体， "系列 2，， 表示癌旁细胞群体；
图 10为本发明 53个癌细胞与 8个正常细胞分类关系示意图；
图 11为本发明癌细胞与正常细胞聚类示意图， 横坐标表示第一主成分向量， 纵坐 标表示第二主成分向量； 图 12为本发明单细胞分类装置示意图；
图 13为本发明单细胞分类装置中筛选模块示意图。 发明详细描述
下面详细描述本发明的实施例，所述实施例的示例在附图中示出， 其中自始至终相 同或类似的标号表示相同或类似的元件或具有相同或类似功能的元件。下面通过参考附 图描述的实施例是示例性的， 仅用于解释本发明， 而不能理解为对本发明的限制。 本发明采用新一代测序技术（NGS ), 通过生物信息学方法， 对单细胞基因组进行分析 和研究， 筛选和收集细胞亚群（或孩 ί粒）进行深入的后续研究。 以更高效、 方便地应用于 临床诊断和治疗 （如产前诊断、 胚胎植入前遗传诊断、 个体化治疗、 多点图谱制作、 精子 和卵子的分型、 遗传病诊断、 肿瘤 （如淋巴瘤、 白血病）分型等）、 医学研究 （如自闭症、 神经系统疾病和自体免疫性疾病的研究、 基因组变异率研究、 干细胞研究、 耐药性研究、 药物靶点基因的筛选、 转录组检测、 细胞模型研究、 种群鉴定等）、 考古学研究、 法医学检 测中。
本发明中涉及的单细胞样本包括核酸（基因组 DNA或 RNA,如非编码 RNA、 mRNA ); 单细胞来自生物体， 采用常规方法制备。 特别的， DNA或 RNA可由细菌、 原生动物、 真 菌、 病毒及高等生物 /高等动植物， 如哺乳动物， 特别是人类的单细胞提取或扩增得到。 单 细胞可经过体外培养， 或者从临床样本 （包括血浆、 血清、 脊液、 骨髓、 淋巴液、 腹水、 胸腔积液、 口腔液体、 皮肤组织， 呼吸道、 消化道、 生殖道、 泌尿道， 眼泪、 唾液、 血细 胞、 干细胞、 肿瘤） 中直接分离得到， 胎儿细胞可来自胚胎 （如一个或几个拟胚 /胚胎）或 母体血液， 可来自活体或者死亡生物体。 样本包括单细胞悬液、 石蜡包埋组织切片、 穿刺 活检组织。
样本可反映细胞特定状态， 如细胞增殖、 细胞分化、 细胞凋亡 /死亡、 疾病状态、 外界 刺激状态、 发育阶段。
特别地， 单细胞样本来自于哺乳动物， 包括植入前胚胎、 干细胞、 疑似癌细胞、 病源 生物， 甚至来自犯罪现场。 例如分析人类卵裂球细胞（八细胞阶段胚胎及以后）， 可检测胎 儿基因组是否出现遗传缺陷。
下面结合图 3 , 对本发明单细胞分类方法的具体实现过程进行详细说明。 其中， 图 3示 出了从（7 )开始的流程。
( 1 )分离单细胞： 采用物理机械、 化学、 生物的方法， 如 流控、 口吸分离、 梯度稀 释、 低熔点琼脂糖固定等方法， 分离得到包含完整基因组的单个细胞。
( 2 ) 裂解细胞： 对分离得到的单个细胞， 采用去污剂法、 煮沸法、 碱变性法、 溶菌酶 法、 有机溶剂法等方法， 裂解细胞核， 得到完整的细胞基因组 DNA。
( 3 )单细胞全基因组扩增 （ WGA ):
目前全基因组扩增有 2种策略： 即基于 PCR的扩增， 如 DOP-PCR、 PEP-PCR、 T-PC , 和线性 DNA扩增， 如 OmniPlexWGA、 多重置换扩增（MDA )。 进行单细胞全基因组扩增， 以达到新一代测序技术所要求的 DN A起始量。
( 4 )全基因组扩增产物定量： 采用凝胶电泳检测、 Agilent 2100 Bioanalyzer检测、 Quant-iT™ dsDNA B Ass^ Kit检测等方法对单细胞全基因组的扩增产物进行定量，结果显 示为无降解、 符合新一代测序技术所要求的 DNA起始量的样本才可以继续 DNA文库构建 以及上机测序。
( 5 )全基因组扩增产物检测： 采用 STR检测、 Housekeeping Gene检测等方法， 对单 细胞全基因组的扩增产物进行检测， 结果显示为扩增产物在相应物种的染色体上均匀分布 的样本才可以继续 DNA文库构建以及上机测序。
( 6 ) DNA 文库构建及上机测序： 采用常规的全基因组 DNA 文库构建或外显子组
( Exo me )序列捕获技术进行 DNA文库构建， 质检合格后 , 采用新一代测序技术进行单细 胞基因组测序， 如 Illumina HiSeq 2000测序系统、 Illumina Genome Analyzer II测序系统、 AB SOLiD™ 4.0测序系统、 Roche GS FLX Titanium System等。
( 7 ) Reads定位
每个单细胞样本经测序得到的 Reads结果用短序列比对软件 (如 SOAPaligner、 BWA、
Bowtie )与参考基因组序列 （如人类基因组 HG18、 HG19 )进行比对， 并根据具体数据设 置最优参数， 以对 Reads进行准确定位。
( 8 )基本数据统计
根据上述比对结果， 计算出每个样本 （单细胞或组织）相对于参考基因组序列的测序 深度和覆盖度等结果。
测序深度 ( sequencing depth )指一个基因组被测到的平均深度， 计算方法为测序的总 碱基数除以基因组的大小。
覆盖度（sequencing coverage )指一个基因组中被测到的部分所占的大概比例， 计算方 法为基因组被覆盖到的位点除以基因组的有效长度。
测序深度和覆盖度用来评价用于生物信息分析的数据量是否足够， 测序是否均一。
( 9 )数据过滤
当存在多对重复双末端 reads时， 在多对 reads的序列完全一致时随机选择一对 reads , 将所述多对 reads中其他重复双末端 reads去除； 和 /或去除未唯一比对上参考基因组序列的 reads
根据数据特点挑选出每个 DNA文库中重复双末端（ p air-end ) reads , 比如， 由于 PCR 扩增次数过多引起的重复双末端 reads, 当然， 不限于 PCR扩增， 也可以是其他扩增方式导 致具有多对重复的双末端 reads。
当存在多对重复双末端 reads, 且所述重复双末端 reads 的序列完全一致时， 从中随机 选择一对， 将其他重复双末端 reads去除。
如图 1所示， 三对 reads A、 B、 C序列完全一致并且比对到基因组上的起止位置也完全 一致， 起止位置完全一致即具有重复双末端。 在这种情况下只随机保留其中一对 reads, 将 其他重复的 reads去除。
为了确认数据的准确性， 还可以去除未唯一比对上参考基因组序列的 reads。 以对人类 基因组进行外显子组测序为例， 当然， 并不限于此， 比如对其他哺乳动物等进行测序， 测 序方式也不限于外显子组测序， 比如全基因组测序等方式。 考虑到人类外显子区域在基因 组上并不可能有多个拷贝， 即不可能是重复序列， 所以外显子组测序得到的 reads应大都唯 一比对上人类参考基因组。 为排除错误比对造成的影响， 只挑选出唯一比对上参考基因组 的 reads进行分析（即 hit数为 1的 reads ), 这样在极大程度上减少了错误比对造成的影响。
如图 2所示， Reads D比对上参考基因组的多个位置， 而 Reads E只比对上唯——个位 置， 由于外显子组在基因组上并不是重复区域， 因此直接去除 Reads D。
( 10 )个体基因型 （ notype ) 的判定
我们充分考虑并利用参考基因组上的已有信息， 用基因型判定软件（如 SOAPsnp、 SAMtools Pileup pileup ), 判定每个细胞样本目标区域的可能的基因型， 以得到每个细胞 样本的基因型的可能性文件。
本发明中所测的是外显子区域的数据， 在该实施例中， 目标区域是外显子所在的区域。 一般会指明所需测序和生物信息分析的特定区域如：
chrl 20038 20358
chrl 58832 59992
chrl 357410 358570 ( 11 ) SNP数据集
由于每个细胞基因组存在一些低深度的区域， 本发明综合所有细胞基因型的可能性文 件， 利用最大似然法（ maximum likelihood approach )整合所有细胞的数据， 产生包含所有 样本每个位点的伪基因组（ Pseudo-genome )。选择概率最大的基因型作为每个细胞的一致基 因型， 通过基因型和测序深度等信息检测出高质量的 SNPs。 得到样本的一致性序列后， 结 果保存为群体 SNPs格式的 SNP数据集。
( 12 )群体 SNP的基因型
根据参考基因组 SNP数据集位置信息， 从基因型的可能性文件中提取每个细胞相应位 点的基因型， 得到群体 SNP 在各细胞相应位点的基因型文件。 位点指 SNP 所在的位置 ( osition )。
( 13 )挑选与细胞突变相关的 SNP位点
本发明主要在于寻找各细胞间的差异位点， 因此须挑选出与细胞突变有关的位点。 首先， 去除非群体间的 SNP位点。
非群体间的 SNP位点的定义： 所有个体碱基型都一样， 且相对参考序列都是 SNP, 那 么该位点是非群体间的 SNP位点。
例如， 参考序列是 A, 所有个体在该位置都是杂合的碱基型 R, 该位点即是非群体间 的 SNP位点。 如 chrl 319660
其次，还可以去除存在杂合丢失的位点。 由于单细胞在 WGA扩增时， 存在只扩增出一 对染色体中的一条的情况， 导致等位基因丢失（Allele Dropout ), 所以检测出来的每个细胞 在某些位点存在杂合丢失的现象。 排除这一类位点的干扰。
最后， 去除已经公布的 SNP位点， 比如， 以人类为例， 去除正常人的 SNP位点， 即去 除人类基因组 HG 18的 dbSNP、 炎黄一号的 SNP和千人基因组的 SNP。
上述的三个操作顺序不分先后， 在执行完这三个操作之后， 得到的 SNP位点是细胞突 变 SNP位点。
( 14 )群体结构分析
根据细胞群体突变 SNP位点的基因型文件，分别运用生物信息群体分析中常用的方法， 如邻接 N-J 方法构树、 MEGA软件、 主成分分析 ( Princi al Com onents Analysis , PCA X 群体结构 （ Group Structure )等对细胞进行分类。 在进行细胞分类时， 可以采用以上方法中 的至少一种。 作为本发明的一个实施例， 采用以上全部方法， 当各种方法的分类结果一致 时， 确认为最终的细胞分类结果。
14-1、 邻接 N-J 方法构树
由于不同类别的细胞受到选择的程度不同， 即单碱基的突变率不同， 因此进化中的类 别差异也在 SNP上表现出来。 两个细胞的差异程度可由 SNP数据计算得到。 经理论推算， 得到如下公式：
上式中 D 表示细胞 和细胞 '差异距离，其中"是 SNP的总数目，而 ^^表示细胞 和 细胞 在位点 的差异程度， 其中定义
基因型完全相同， 如位置 ， 细胞 : ， ':
基因型完全不同， 如位置 ， 细胞 ， ': C
Figure imgf000012_0001
基因型部分不同， 如位置 ， 细胞 :^，j': M 由于人类基因组是二倍体， A表示等位基因中的两个位点都是 A, M是杂合位点， 即 是 A和 C组合。 基于上述步骤（ 13 )获得的细胞群体突变 SNP位点的基因分型文件， 统计 所有单细胞样本两两对比的差别， 得到一个目标区域两两对比的差异矩阵， 将上述差异矩 P车应用 l) Fneighbor程序 ( http：〃 emboss.bioinformatics.nl/cgi-bin/ emboss help/fneighbor ), 该 程序通过邻接（N-J ) 的方法可得到其进化树。
14-2、 MEGA软件
利用 MEGA软件 ( http ://www.megasoftware.net ), 将每个细胞所有 SNP位点的基因型 文件组成序列， 一个细胞对应一个序列， 作为 MEGA的输入文件， MEGA根据各细胞序列 上的差异， 其中该软件有三种方法 ( M aximum likelihood、 Least Squares 和 M aximum p ars imony )构成关系树。
14-3、 PC A主成分分析
在统计学中， 主成分分析 ( Princi al Com onents Analysis, PCA )是一种筒化数据集的 技术， 是一个线性变换。 这个变换把数据变换到一个新的坐标系统中， 使得任何数据投影 的第一大变量数在第一个坐标（称为第一主成分）上、 第二大变量数在第二坐标（第二主 成分）上， 依次类推。 主成分分析经常用于减少数据集的维数， 同时保留对数据集贡献最 大的特征变量。 通过保留低阶主成分， 忽略高阶主成分来实现的。 这是由于低阶成分往往 能够保留数据集中最重要的方面。
根据参考文献 A tutorial on Princi al Components Analysis. Lindsay I Smith, 2002-02和真 实的 SNP数据特点， 首先将 SNP数据转换成数字矩阵（与参考序列一致为 0 , 相反为 2 , 筒并碱基为 1 ) 并作均一化。 然后通过上述介绍的方法构建线性向量方程。
= ai0 +anxt + ai2x + ....+ amx 其中 从 1到 21表示第 i个样本。 应用 R语言软件包强大的解方程能力， 解得矩阵"， 根据各细胞的数据特点提取前四 个主成分向量， 以向量作为坐标轴展示各细胞聚类情况。
14-4、 Structure群体结构
Structure软件 ( htt ://prit ch.bsd.uchicago.edu/software/structure2 1. html )基于 SNP位点 的基因分型数据， 推断是否存在不同群体并判断每个个体所归属的群体。 根据软件说明， 将群体 SNP的基因型文件转换格式， 作为 Structure输入文件并在混合模型中采用高达 5万 次模拟， 假设多个群体存在时， 计算每个细胞归属各类群体的概率。
经过以上的方法流程， 实现了对单细胞的分类。 在分类的基础上， 还可以进一步筛选 单细胞， 其流程如下：
( 15 )群体结构分析结果
根据上述群体结构分析结果， 实现了对单细胞的分类， 提取每个细胞样本信息， 剔除 存在异议的细胞， 比如分类不明确或明显离群样本。
( 16 ) 筛选相关基因
根据细胞群体的 SNP , 通过一系列统计和检验进行这些群体在基因组的比较， 找出有 明显差异的区域或者基因， 即可筛选到相关系数较高的基因。
以人类基因组为例， 具体做法如下：
从人类基因组数据库中下载 HG18 对应的注释文件， 得到目前已预测的人类基因组中 30000多个基因的起止位置。
根据分类结果获得已分类群体， 计算每类群体中每个基因所有 SNP位点的统计量， 并 累加统计量。 这里所说的每个基因是指已预测的基因组中的基因。 其中主要采用的计算统计量 的公式如下， 是衡量一个群体多态性高低的指标， a'b 指某群体中两种碱基的样本数， 公式可以为：
a * b
π =—； ~
C2 h 还可以对获得的统计量作差异检验， 获得检验值。 所采用的检验值至少为如下之一： Lod、 Fst、 Pbs。 作为本发明的一个实施例， 可以采用以上三种检验值， 当以上三种检验 值一致时， 作为最终的检验值结果。
将这 30000多个基因按统计量和 /或检验值进行排序， 挑选出统计量和 /或检验值最高的 基因。 也就是说， 可以根据统计量进行排序、 可以根据检验值进行排序、 也可以根据统计 量和检验值进行排序。 作为本发明的一个实施例， 可以采用最后一种方法， 当根据统计量 得到的排序结果与根据检验值得到的排序结果一致时， 作为最终筛选出的基因。
( 17 )基因功能分析
查看筛选出的基因的功能， 分别进行功能分析。 判断这些基因是否在某些通路中受到 影响， 从而与某些疾病的发生机制相关。
如图 12所示， 为本发明单细胞分类装置示意图。 该装置包括：
数据过滤模块，将每个单细胞样本经测序得到的 reads结果与参考基因组序列进行比对， 并将比对结果进行数据过滤；
基因型确定模块， 与数据过滤模块耦合， 根据过滤后的数据确定每个单细胞样本的一 致基因型， 并将所有单细胞样本的一致基因型保存为 SNP数据集；
基因型文件提取模块， 与基因型确定模块耦合， 从已保存的 SNP数据集提取与参考基 因组 SNP数据集位置对应的位点的基因型文件；
分类模块， 与基因型文件提取模块耦合， 挑选细胞突变 SNP位点， 根据细胞群体突变 SNP的基因型文件， 对细胞进行分类， 采用的分类方法至少包括如下之一： 邻接 N-J 方法 构树、 MEGA软件、 PCA主成分分析、 Structure群体结构。
在另一实施例中， 还如图 12所示， 单细胞分类装置还包括：
筛选模块， 与分类模块耦合， 用于得到已预测基因组中基因的起止位置； 根据分类结 果获得已分类群体， 计算每类群体中每个基因所有 SNP位点的统计量， 并累加统计量； 对 获得的统计量作差异检验， 获得检验值； 将已预测基因按统计量或检验值进行排序， 筛选 出统计量或检验值最高的基因。
筛选模块可以进一步包括以下单元， 如图 13所示， 其中：
获取单元， 得到已预测基因组中基因的起止位置；
计算单元， 根据细胞分类结果获得已分类群体， 计算每类群体中每个基因所有 SNP位 点的统计量， 并累加统计量； 对获得的统计量作差异检验， 获得检验值；
排序单元， 耦合于获取单元和计算单元， 将已预测基因按统计量或检验值进行排序， 筛选出统计量或检验值最高的基因。
本发明单细胞分类装置中各个模块执行的具体操作体现在上述方法流程中， 也可以根 据以下实施例看出各个模块的具体操作。 下面将结合实施例对本发明的方案进行解释。 本领域技术人员将会理解， 下面的实施 例仅用于说明本发明， 而不应视为限定本发明的范围。 实施例中未注明具体技术或条件的， 按照本领域内的文献所描述的技术或条件 (例如参考 J.萨姆布鲁克等著， 黄培堂等译的《分 子克隆实验指南》， 第三版， 科学出版社）或者按照产品说明书进行。 所用试剂或仪器未注 明生产厂商者， 均为可以通过市购获得的常规产品， 例如可以采购自 Illumina公司。
实施例 1: 腎癌单细胞的分类
1-1、 Reads定位
每个单细胞样本经测序得到 的 reads 结 果用 SOAPaligner 比对软件 ( htt ://soa .genomics.org. cn/soa aligner.html )与参考基因组序列（人类基因组 HG18 )进行 比对， 由于人类 SNP在千份之二并且 Reads的读长为 lOObp , 所以 SOAP比对时设置每条 Reads最多有 3个错误比对（ M ismacth ), 不容 Gap的参数， 以确保可比对上的 Reads的位 置准确。
1-2、 基本数据统计
根据上述比对结果， 计算出每个样本 （单细胞或组织）相对于参考基因组序列的测序 深度和覆盖度等结果， 经统计得到全基因组测序且 Mean Depth为 3x附近时， 由于 PCR扩 增存在一定的偏向性（ Bias ), 所以样本的覆盖度在 55-90%之间波动较大。
表 1 肾癌单细胞样本全基因组测序的覆盖度和深度数据
Figure imgf000015_0001
其中 RC-1至 RC-15表示肾癌单细胞， 共有 15个单细胞样本； RN-1至 RN-5表示癌旁 单细胞； RN-T表示正常组织直接提取 DNA测序， 以作为对照用于数据分析和评估。 主要 以癌旁单细胞作为对照样本。 也存在以癌旁单细胞以及正常组织都同时作为对照样本， 如 在去除杂合丢失的位点时就用到上述两个对照样本。
在外显子组测序中增加测序深度 , 当外显子目标区域的平均深度 ( M ean Depth )为 30x 附近时， 目标区域覆盖度达到 80-96%。 从统计学意义上看， 一个位点有四条 reads 支持则 可以判断该位点碱基的准确性达到 99%, 而统计得到测序深度为 4的位点占外显子区域的 比例达到 60-90% , 说明外显子组测序的数据比全基因组测序得到的数据更好。 表 2 肾癌单细胞样本外显子组测序的覆盖度和深度数据
Figure imgf000016_0001
将上述 2 个表格进行对比， 可得到全基因组测序的深度低而无法进行后续分析， 而外 显子组测序的深度高。 另外考虑到测序成本问题， 因此下面主要基于外显子组测序得到的 数据进行分析。
1-3、 数据过滤
根据数据特点挑选出每个 DNA文库中由于 PCR扩增次数过多引起的重复双末端 reads , 在多对重复双末端 reads的序列完全一致时随机挑选一对 reads , 去除其他 reads。
比如， 图 1中三对 reads A、 B、 C序列完全一致并且比对到基因组上的起止位置也完全 一致。 在这种情况下只随机保留其中一对 reads。
为了确认数据的准确性， 考虑到人类外显子区域在基因组上并不可能有多个拷贝， 即 不可能是重复序列， 所以外显子组测序得到的 reads应大都唯一比对上人类参考基因组。 为 排除错误比对造成的影响，只挑选出唯一比对上参考基因组的 reads进行分析（即 hit数为 1 的 reads ), 这样在极大程度上减少了错误比对造成的影响。
如图 2 , Reads D 比对上参考基因组的多个位置， 而 Reads E只比对上唯——个位置， 由于外显子在基因组上并不是重复区域， 因此直接去除 Reads D。 1-4, 个体基因型的判定
我们充分考虑并利用人类基因组（该实施例中为参考基因组）上的已有信息， 从 NCBI 网站下载人类基因组（HG18 )对应的 dbsrp, 作为 SOAPsrp的先验概率， 并依据目前确定 的研究结果，设置杂合位点 SNP的先验概率为 0.0010,纯合位点 SNP的先验概率为 0.0005。
在设置以上参数后， 将步骤 1-3经过滤的数据输入到 SOAPsnp软件中， 由 SOAPsrp软 件将过滤后数据与参考基因组比对， 得到比对结果为 CNS文件。
1-5、 SNP数据集
由于每个细胞基因组存在一些低深度的区域， 本发明综合所有细胞基因型的可能性的 文件， 利用最大似然法整合所有细胞的数据， 产生包含所有样本每个位点的伪基因组 ( Pseudo-genome )„选择概率最大的基因型作为每个细胞的一致基因型，通过基因型和深度 等信息检测出高质量的 SNPs。 得到样本的一致性序列后， 结果保存为群体 SNPs格式。
1-6、 SNP的基因型
根据参考基因组 SNP数据集位置信息，从 CNS文件中提取每个细胞基因型的相应位点， 得到群体 SNP在各细胞相应位点的基因型文件。 文件格式如表 3所示。
"SNP位置"表示该 SNP位点在染色体上的位置， "碱基型"对应某一个细胞基因组在 这一位点的碱基型， 深度为 0的位点用 "-"表示（即数据丢失的位点）。 "样本 ID"对应 21 个单细胞或组织 DNA样本。
表 3 群体 SNP在各细胞相应位点的基因型文件格式示意
Figure imgf000017_0001
该 21个单细胞或组织 DNA样本相对于人类基因组 HG18在目标区域的群体 SNP数目 为 93957个。 其中杂合位点的组合由以下字母表示：
"M"代表 "A和 C" , " "代表 "A和 G" , "W"代表 "A和 T" , "Y"代表 "C和 T" , "S" 代表 "C和 G" , "K" 代表 "G和 T"。
1-7、 挑选与细胞突变相关的 SNP位点
本发明主要在于寻找各细胞间的差异位点， 因此须挑选出与细胞突变有关的位点。
表 4 非群体间 SNP位点示意
Figure imgf000017_0002
chrl 10402265 R R R R R R R R - R R R R R R R R R R R R chrl 11001664 R R R R R R R R R R R R R R R R R R R R R chrl 12775804 W W W W W W W W W W W W W W W W W W W W W chrl 12775818 Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y Y 首先， 去除非群体间的 SNP位点， 这些位点如表 4所示。 在同一位点所有单细胞基因 组的碱基型都一致， 即在该位点 21个细胞构成的群体都是一致的。 经过计算， 共有 504个 这样的位点， 去除之后剩下 93453个 SNP位点。
其次， 由于单细胞在 WGA扩增时， 存在只扩增出一对染色体中的一条的情况， 导致等 位基因丢失（Allele Dropout ), 即原本是杂合位点， 在测序时只测到其中一个碱基型， 所以 检测 1出来的每个细胞在某些位点存在杂合丢失的现象， 如表 3所示的位点：
chrl 9 10079332 G
chrl 9 10079408 R R R R R R R R G R R R R R R R R R R R R
即第 9个单细胞样本 RC-9在该位点很可能原本是杂合位点，但由于只扩增到一个碱基 型， 所以判断为纯合位点。
为了排除这一类位点的干扰， 考虑到同一位点同时出现几个样本都发生杂合丢失的概 率极低 , 因此本发明采取策略如下：
第一， 正常组织直接提取 DNA测序的 RN-T (即最后一列） 必须为杂合位点。 因为这 个样本不经过 WGA扩增， 不可能出现杂合丢失。
第二， 其它 20个单细胞样本中杂合位点和数据丢失的样本数大于等于 18。 即同一位点 最多只允许两个单细胞样本发生杂合丢失。 因为三个或三个以上样本在同一位点同时发生 杂合丢失的概率极低。
上述两个条件都要满足， 即正常组织直接提取 DNA测序的最后一列为杂合位点， 其它 20个单细胞样本中杂合位点和数据丢失的样本数大于等于 18。 这样的位点才除去。
经过计算， 共有 3975个这样的位点， 经过这一步过滤得到 SNP数目为 89478个。 最后， 为了得到与肾癌单细胞突变有关的位点， 须去除已经公布的正常人的 SNP位点， 即去除人类基因组 HG18的 dbSNP、 炎黄一号的 SNP和千人基因组的 SNP, 得到与各细胞 突变相关的 SNP位点 50524个。
1-8、 群体结构分析
根据细胞群体突变 SNP位点的基因型文件， 分别运用生物信息群体分析中常用的方法 对细胞进行分类。 分类是由进化树的分枝和聚类情况来确定的。 如图 4, RC和 RN在进化 树中明显聚集为分离的两部分， 因此将其分为两类。
1-9.1 邻接 N-J 方法构树
如图 4所示， 根据进化树可以对细胞进行分类。
1-9.2 MEGA软件
图 5是最大似然法构成的关系树， 根据关系树对细胞进行分类。
1-9.3 PCA主成分分析 如图 6、 图 7和图 8所示的肾癌外显子组测序 PCA结果， 根据聚类情况对细胞进行分 类。
1-9.4 Structure群体结构
如图 9, 横坐标表示样本名称， 纵坐标表示每个样本所属群体的概率， 根据概率对单细 胞进行分类。 由图 9所示， 20个单细胞大致可以分为两个群体。 肾癌外显子组测序 Structure 结果如图 9所示。
1-10、 群体结构分析结果
根据上述群体结构分析结果， 提取每个细胞样本信息， 剔除存在异议的细胞（分类不 明确或明显离群样本）。 从上述各种群体结构分析结果来看， 取样正常， 分类合理。 这 20 个单细胞样本大致可以分为 2个群体， 即癌细胞群体（ 15个 RC )和癌旁细胞群体（5个 RN ), 其中 RC-1、 RC-6、 RC-7为癌细胞的一个亚群。
细胞样本信息指进行分析的单细胞中哪些是癌细胞、哪些是癌旁细胞（由取样时确定）， 细胞样本信息仅作为参考， 需要与聚类结果一起分析。 若取样时的细胞样本信息认为是癌 细胞和癌旁细胞， 而聚类时正好分成两个群体聚集， 说明结果互相对应； 若取样时的细胞 样本信息与聚类结果不一致， 以聚类结果为准。
在癌细胞聚类中另行聚类在一起， 由此确认 RC-1、 RC-6、 RC-7为癌细胞的亚群。
1-11、 筛选与肾癌相关的基因
根据上述 RC和 RN两个细胞群体在外显子区域 SNP,通过一系列统计和检验进行这两 个群体在外显子区域的比较， 找出有明显差异的区域或者基因， 即可筛选到与该例肾癌相 关系数较高的基因。 具体做法如下：
1-11.1 从人类基因组数据库中下载 HG18 对应的注释文件， 得到目前已预测的人类基 因组中 30000多个基因的起止位置。
1-11.2根据分类结果， 得到 RC和 RN两个群体， 计算每类群体中每个基因所有 SNP 位点的统计量， 并累加所述统计量。
其中主要采用的计算统计量 的公式如下， 是衡量一个群体多态性高低的指标， a,b 指某群体中两种碱基的样本数， 公式可以为：
a * b
π =― ~
C2 h 如在上述 15个 RC样本中， 共有 30条染色体， 对于下述两个位点： 位点 1 只有 1条 染色体为 C, 而另外 29条染色体为 T ( a=l , b=29 ); 位点 2有 15条染色体为 C, 而另外 15条染色体为 T (a=15 , 15)。代入公式，得到位点 1的 值是 0.06,位点 2的 值是 0.517, 则这 2个位点的多态性具有明显差别。
在统计一个基因的多态性时， 累加该基因所有位点的 值， 由于非 SNP位点的 值为
0 ( " = 0或 6 = 0时， ^=0 ), 即对某群体累加该基因所有 SNP位点的 值。 1-11.3 将这 30000多个基因按统计量或检验值进行排序， 挑选出统计量或检验值最高 的基因， 查看这些基因的功能。
所采用的检验值至少为如下之一： 检验值为^^ 、 Fst、 Pbs , 该实施例采用以上三种 检验值。 下面具体说明各个检验值以及计算过程。
分别对 RC和 RN这两个群体代入数据， 可以得到 ^ ^和 5^!。 由于需要对这两个群 体的差异进行比较， 因此定义^^如下：
Figure imgf000020_0001
若 Αϋ和 相差很小， 即该基因在这两个群体中没有太大差异， Lod( 可 明显偏离于 0 , 则可初步认为该基因是造成这两个群体分化的重要基因。 如上所述， 分别统计 HG18中 30000多个基因的值 £c¾ g^ , 然后从大到小进行排序， 筛选出排名较前的基因。
FsT ( Fixation index)主要用来评价群体间的基因组距离和种群的差异， 量种群间 分化程度的一个指标， 由 Sewall Wright在 1922年应用 F-检验的一种特殊情况发展而来。
FsT的零假设是在群体没有分化时， 多态性位点在群内和群间的次等位碱基的频率差别 是不显著的。 计算 F^的方法很多， 虽然具体计算方法不同， 但基本理论是一致的， 即由 Hudson ( 1992 )给出的定义：
Z _ ^Between― ^Within
^Betwe en 其中11^1表示从两个群体（Between )中分别抽取一个样本， 组成一对， 计算这对样 本 SNP基因型的差异， 如此可以计算所有成对样本 SNP基因型的差异， 最后求平均值即为
^Between
Π ί η表示从一个群体（Within )中分别抽取 2个样本， 组成一对， 计算这对样本 SNP 基因型的差异， 如此可以计算所有成对样本 SNP基因型的差异， 最后求平均值即为 nwi in。 如果有两个群体， 是两个群体分别先计算11 然后累加。
结合目前已有 SNP集的数据结构， 基于上述原理， 推导公式如下：
Figure imgf000021_0001
上式中 X是 SNP 位点 i在群体 j中的次等位碱基 (第二碱基）的频率； 而 是 SNP 位 点 在群体 · 中染色体上的物理位置； 则是群体 用于分析的 SNP位点个数的总和。 其中变量 取 RC和 RN , 变量 以最后判定的 SNP位置代入。 以基因为单位计算每个 基因的 值， 然后对 HG18中 30000多个基因的 值进行排序， 筛选出排名较前的 基因。
由于数据丢失（Missing data ) 的情况下， 对 SNP位点频率估计不精确， 使得 Fsr不能 灵敏反映出数据原本的属性。 根据参考文献 Sequencing of 50 Human Exomes Reveals Adaptation to High Altitude. Science, 2 July 2010, 329, 75-78 )采用的方法， 对 取 log, 并引 入第三个群体（本实施例引入千人基因组中的部分数据， 北京人基因组数据记为 B ), 定义
Ί Pbs , 公式如下：
T = - \og(l - Fst) 即三个群体两两对比的^ ^如下：
Figure imgf000021_0002
TRC-B = - ^g(\ - FstRC_B )
此时^ w的公式如下： pbs― TRC RN + TRC—B ― TRN—B
2 以基因为计算单位， 分别统计 HG18中 30000多个基因的 值， 然后排序筛选出 排名较前的基因。
1-12、 基因功能分析 根据以上三种检验值中至少一个，该实施例根据 ^^以及尸&筛选得到重要基因， 分别进行功能分析。 判断这些基因是否在某些通路中受到影响， 从而与肾癌的发生机制相 关。 实施例 2: 白血病单细胞的分类和筛选
2-1、 Reads定位
对每个癌症单细胞进行 30x深度的外显子组测序， 得到的 reads结果用 SOAPaligner2.0 比对软件与参考基因组序列 （人类基因组 HG18 )进行比对。 由于人的 SNP在千份之二并 且 Reads的读长约为 lOObp , 所以 SOAP比对时我们设置每条 Reads最多有 2个错误比对 ( mismacth ), 不允许出现 Gap , 以确保比对到参考基因组上 Reads的准确性。
2-2、 基本数据统计
共进行 53个癌细胞、 8个口腔上皮细胞（正常细胞）测序。 表 5为各细胞样本外显子 组测序的覆盖度和深度数值信息。
表 5 各细胞样本外显子组测序的覆盖度和深度
Figure imgf000022_0001
ET-22 0.85 24.00 ET-21 0.63 16.00
ET-6 0.85 17.00 ET-27 0.62 17.00
ET-87 0.85 34.00 ET-15 0.60 16.00
ET-16 0.84 23.00 NC-30 0.46 22.00
ET-4 0.84 15.00 NC-7 0.32 6.01
ET-43 0.84 25.00 NC-17 0.29 15.00
ET-5 0.84 17.00 NC-29 0.25 8.62
ET-25 0.83 20.00 NC-5 0.24 4.65
ET-94 0.83 40.00 NC-28 0.21 4.06
ET-3 0.81 23.00 NC-14 0.21 5.40
ET-91 0.81 31.00 NC-8 0.21 5.85
ET-29 0.80 18.00
2-3、 数据过滤
与实施例 1相同
2-4、 个体基因型的判定
与实施例 1相同
2-5、 SNP数据集
在进行 SNP数据集的确定时， 考虑到血癌细胞数目较多， 每个单细胞基因组外显子的 覆盖率不是很高， 且是基于每个个体进行 SNPs的确定， 所以我们选用较为严格的标准对得 到的数据进行筛选。
标准如下：
Soap snp软件中每个位点的一致性基因型的质量值不小于 20 , 秩检验 p值不小于 1 %; 对于杂合变异的 SNPs: 位点基因型与参考基因组不同， 最好碱基型测序质量值不小于 20、 测序深度不小于 6, 次好碱基型测序质量值不小于 20, 测序深度不小于 2, 且两基因型测序 深度的比值在 0.2-5范围内。
质量值越大，基因分型越正确，一般大于 20时，错误率在万分之一以下， 可忽略不计。 利用如上标准筛到可靠的 SNPs以后， 根据参考基因组 SNP数据集位置信息确定位点， 提取每个细胞每个位点的基因分型数据生成基因型文件。 文件格式如表 3所示。
2-6、 群体结构分析
根据细胞群体突变 SNP的基因型文件， 我们分别运用多种生物信息群体分析中常用的 方法对各细胞分类。
2-6.1、 邻接 N-J 方法构成聚类构树
如图 10所示为本发明 53个癌细胞与 8个正常细胞分类关系示意图， 其中， ET-T 1表示 癌组织， NC-T 1表示正常组织。
2-6.2、 PCA主成分分析 如图 11所示为本发明癌细胞与正常细胞聚类示意图， 其中， LC表示癌细胞， LN表示 正常细胞。
根据上面的群体分析结果， 提取细胞样本信息， 剔除存在异议的细胞（分类不明确， 或离群样本）。 从上面的群体结构说明取品正常， 分类合理。
2-6.3、 亚群的分类
依据进化树的形状或者情况， 在所有 53个癌细胞中可以清晰地分出 4类亚群， 说明在 癌细胞中存在真实的差异。 利用单细胞分析方法可以将同一癌症组织中不同的细胞亚群分 类出来。
2-7、 挑选高可信度的体细胞突变
从基因型文件中筛高可信度的体细胞突变位点， 标准如下：
正常细胞拥有一致的纯合基因型， 癌细胞存在两个或两个以上的杂合突变或纯合突变， 且不能出现第三种纯合基因型以及和两种纯合基因型不一致的杂合基因型。 比如正常细胞 基因型为 A, 或突变类型是 A->C, 则癌细胞中只能出现三种基因型， 即 、 C和 M, 且 C 和 M 的数量不小于 2。 这种位点我们称为高可信度体细胞位点 （High confident somatic mutation, HCSM )。由于我们利用的是外显子组测序技术，所以过滤掉非外显子区域的位点， 共得到 2296个 HCSMs, 其中同义位点 879个， 非同义位点 1417 (包含错义突变和截断突 变位点）， 非同义 /同义突变比值为 1.61 , 表 6所示。
表 6 高可信度的体细胞突变位点
Figure imgf000024_0001
2-8、 基因功能和途径分析
是在进行细胞分类和筛选之后可进行的下游分析， 将基因突变位点的位置信息以及每 个基因中存在的非同义突变位点数目作为基因功能富集的标准，利用 Webgestalt在线分析工 具 对 突 变 影 响 基 因 的 功 能 以 及 途 径 进 行 研 究 ( http://bioinfo.vanderbilt.edu/webgestalt/option.php ),发现突变主要集中在如下 8类功能的基 因中。
表 7 突变位点基因功能分析结果
Figure imgf000024_0002
molecular GO:0005198 structural molecule activity 0.0021 function
cellular GO:0044430 cytoskeletal part 0.0003 component
cellular GO:0043228 non-membran e-bounded organelle 0.0008 component
cellular GO:0043232 intracellular non-membrane-bounded 0.0008 component organelle
cellular GO:0005856 cytoskeleton 0.0008 component
cellular GO:0044420 extrac ellu lar m atrix p art 0.0044 component
途径分析发现突变基因主要集中在 10种途径中，其中绝大多数与癌症的发生机制相关：
Metabolic athw¾^s (代谢途径 )
ECM -receptor interaction ( ECM受体才目互作用 )
Pathways in cancer (癌症涉及的通路 )
Viral myocarditis (病毒' |·生心月几炎）
Ty e I diabetes mellitus ( I型糖尿病）
MAPK signaling athwi^ ( MAPK信号通路）
Focal adhesion ( |占附斑）
Pantothenate and CoA biosynthesis (泛酸和辅酶 A的生物合成 )
Cell adhesion molecules (CAMs , 细胞黏附分子）
Allograft rejection (移植排斥）
2-9、 突变位点基因功能的预测
我们将外显子区域中非同义突变位点挑出， 利用 SIFT ( http://sift.jcvi.or^ )软件对这些 突变位点对应的基因功能进行预测。 结果分为 4种情况， 如下表 8所示：
表 8 突变位点基因功能预测结果
Figure imgf000025_0001
4种情况是： 对基因功能存在危害； 低可信度危害， 变异对于基因的功能影响不大以及 无法判断。 我们挑选含有功能危害性的突变位点并在上述功能富集和途径中出现的基因进 行后续实验的验证。 工业实用性
本发明的技术方案， 能够有效地应用于细胞分类和变异发生相关基因的筛选。 尽管本发明的具体实施方式已经得到详细的描述， 本领域技术人员将会理解。 根据已 经公开的所有教导， 可以对那些细节进行各种修改和替换， 这些改变均在本发明的保护范 围之内。 本发明的全部范围由所附权利要求及其任何等同物给出。
在本说明书的描述中， 参考术语 "一个实施例"、 "一些实施例"、 "示意性实施例"、 "示 例"、 "具体示例"、 或 "一些示例" 等的描述意指结合该实施例或示例描述的具体特征、 结 构、 材料或者特点包含于本发明的至少一个实施例或示例中。 在本说明书中， 对上述术语 的示意性表述不一定指的是相同的实施例或示例。 而且， 描述的具体特征、 结构、 材料或 者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。

Claims

权利要求书
1、 一种单细胞分类方法， 包括以下步骤：
对来自相同群体的多个单细胞样本的全基因组分别进行测序， 以便获得来自每个单细 胞样本的 reads;
分别将来自每个单细胞样本的 reads与参考基因组序列进行比对， 并对所述 reads进行 数据过滤；
基于经过过滤的 reads , 确定每个单细胞样本的一致基因型， 其中， 所有单细胞样本的 一致基因型构成所述群体的 SNP数据集；
针对所述每个单细胞， 基于所述群体的 SNP数据集， 确定在与参考基因组 SNP数据集 的位置相应的位点， 每个细胞所对应的基因型； 以及
选择与细胞突变相关的 SNP位点， 并且基于所述单细胞在该位点的基因型， 对所述单 细胞进行分类。
2、 根据权利要求 1所述的单细胞分类方法， 其特征在于， 利用第二代或者第三代测序 平台进行所述测序，
其中， 所述数据过滤的标准为：
当存在多对重复双末端 reads时， 在多对 reads的序列完全一致时随机选择一对 reads , 并将所述多对 reads中其他重复双末端 reads去除； 和 /或
去除未唯一比对上所述参考基因组序列的 reads。
3、 根据权利要求 1或 2所述的单细胞分类方法， 其特征在于， 基于经过过滤的 reads , 确定每个单细胞的一致基因型进一步包括：
基于所述经过过滤的 reads , 确定每个单细胞样本在目标区域基因型的可能性； 基于所有单细胞样本在目标区域基因型的可能性， 确定包含所有样本每个位点的伪基 因组； 以及
从所述伪基因组中选择概率最大的基因型作为每个单细胞样本的一致基因型。
4、 根据权利要求 1、 2或 3所述的单细胞分类方法， 其特征在于， 选择与细胞突变相 关的 SNP位点进一步从所述群体 SNP数据集中去除下列的至少之一：
非群体间的 SNP位点、 存在杂合丢失的位点、 以及已经公布的 SNP位点。
5、 根据权利要求 4所述的单细胞分类方法， 其特征在于， 所述多个单细胞样本至少之 一的全基因组在进行测序之前被进行全基因组扩增处理， 其中，
去除存在杂合丢失的位点进一步包括去除符合下列条件的位点：
针对未经全基因组扩增的样本， 测序结果为杂合位点； 并且
在经过全基因组扩增的样本中， 在同一位点， 杂合位点和数据丢失的样本数目大于等 于经过全基因组扩增的样本的数目减去 3。
6、 根据权利要求 1至 5中任一项所述的单细胞分类方法， 针对所述每个单细胞， 基于 所述群体的 SNP数据集， 确定在与参考基因组 SNP数据集的位置相应的位点， 每个细胞所 对应的基因型进一步包括按照下列标准对所述 SNP数据集进行筛选：
每个位点的一致性基因型的质量值不小于 20 , 秩检验 p值不小于 1 %； 以及
对于杂合变异的 SNPs: 最好碱基型测序质量值不小于 20、 测序深度不小于 6, 次好碱 基型测序质量值不小于 20,测序深度不小于 2,且两基因型测序深度的比值在 0.2-5范围内。
7、 根据权利要求 1至 6中任一项所述的单细胞分类方法， 其特征在于， 对细胞进行分 类之后， 还包括以下步骤：
提取每个细胞样本信息， 剔除存在异议的细胞。
8、 根据权利要求 1或 7所述单细胞分类方法， 对所述单细胞进行分类之后， 进一步包 括：
基于分类结果确定已分类群体， 并计算每类群体中每个基因所有 SNP位点的统计量， 任选地对所获得的统计量作差异检验， 获得检验值；
选择统计量或检验值最高的基因或者群体。
9、 一种单细胞分类装置， 其特征在于， 包括：
数据过滤模块， 所述数据过滤模块适于将分别将来自每个单细胞样本的 reads与参考基 因组序列进行比对， 并对所述 reads进行数据过滤， 其中， 所述每个单细胞样本的 reads是 通过对多个单细胞样本的全基因组分别进行测序获得的；
基因型确定模块， 所述基因型确定模块适于基于经过过滤的 reads , 确定每个单细胞样 本的一致基因型， 其中， 所有单细胞样本的一致基因型构成所述群体的 SNP数据集；
基因型文件提取模块， 所述基因型文件提取模块适于针对所述每个单细胞， 基于所述 群体的 SNP数据集， 确定在与参考基因组 SNP数据集的位置相应的位点， 每个细胞所对应 的基因型；
分类模块， 所述分类模块适于基于预先选择的与细胞突变相关的 SNP位点， 并且基于 所述单细胞在该位点的基因型， 对所述单细胞进行分类。
10、 根据权利要求 9 所述单细胞分类装置， 其特征在于， 所述数据过滤模块适于基于 下列标准进行数据过滤：
当存在多对重复双末端 reads时， 在多对 reads的序列完全一致时随机选择一对 reads , 并将所述多对 reads中其他重复双末端 reads去除； 和 /或
去除未唯一比对上所述参考基因组序列的 reads。
11、 根据权利要求 9或 10所述单细胞分类装置， 其特征在于， 所述基因型确定模块适 于通过下列确定所述每个单细胞的一致基因型：
基于所述经过过滤的 reads , 确定每个单细胞样本在目标区域基因型的可能性； 基于所有单细胞样本在目标区域基因型的可能性， 确定包含所有样本每个位点的伪基 因组； 以及
从所述伪基因组中选择概率最大的基因型作为每个单细胞样本的一致基因型。
12、 根据权利要求 9、 10或 11所述的单细胞分类装置， 其特征在于， 分类模块适于从 所述群体 SNP数据集中去除下列的至少之一， 选择与细胞突变相关的 SNP位点： 非群体间的 SNP位点、 存在杂合丢失的位点、 以及已经公布的 SNP位点。
13、 根据权利要求 12所述单细胞分类装置， 所述多个单细胞样本至少之一的全基因组 在进行测序之前被进行全基因组扩增处理， 其中， 所述分类模块适于去除符合下列条件的 位点， 以便去除存在杂合丢失的位点：
针对未经全基因组扩增的样本， 测序结果为杂合位点； 并且
在经过全基因组扩增的样本中， 在同一位点， 杂合位点和数据丢失的样本数目大于等 于经过全基因组扩增的样本的数目减去 3。
14、 根据权利要求 9至 13中任一所述的单细胞分类装置， 其特征在于， 所述基因型文 件提取模块适于按照下列标准对所述 SNP数据集进行筛选：
每个位点的一致性基因型的质量值不小于 20, 秩检验 p值不小于 1%; 以及
对于杂合变异的 SNPs: 最好碱基型测序质量值不小于 20、 测序深度不小于 6, 次好碱 基型测序质量值不小于 20,测序深度不小于 2,且两基因型测序深度的比值在 0.2-5范围内。
15、 根据权利要求 9至 14中任一所述的单细胞分类装置， 其特征在于， 所述分类模块 进一步适于提取每个细胞样本信息， 剔除存在异议的细胞。
16、根据权利要求 9或 15所述的单细胞分类装置， 其特征在于， 进一步包括筛选模块： 基于分类结果确定已分类群体， 并计算每类群体中每个基因所有 SNP位点的统计量， 任选地对所获得的统计量作差异检验， 获得检验值； 以及
选择统计量或检验值最高的基因或者群体。
17、 一种基因筛选方法， 包括以下步骤：
根据权利要求 1-8任一项所述的方法， 对细胞进行分类， 以便获得已分类亚群体， 并计 算每类亚群体中每个基因所有 SNP位点的统计量， 任选地对所获得的统计量作差异检验， 获得检验值；
选择统计量或检验值最高的基因， 作为与细胞突变有关的基因。
18、 一种基因筛选装置， 包括：
细胞分类装置， 所述细胞分类装置为权利要求 9-16任一项所述， 以便对细胞进行分类， 获得已分类亚群体；
计算单元， 所述计算单元适于根据细胞分类结果获得已分类亚群体， 并计算每类亚群 体中每个基因所有 SNP位点的统计量，任选地对所获得的统计量作差异检验， 获得检验值； 排序单元， 所述排序单元将所有基因按统计量或检验值进行排序， 筛选出统计量或检 验值最高的基因， 作为与细胞突变有关的基因。
PCT/CN2012/080578 2011-08-25 2012-08-24 单细胞分类方法、基因筛选方法及其装置 WO2013026411A1 (zh)

Priority Applications (2)

Application Number Priority Date Filing Date Title
US14/239,650 US20140206006A1 (en) 2011-08-25 2012-08-24 Single cell classification method, gene screening method and device thereof
EP12825673.2A EP2749655B2 (en) 2011-08-25 2012-08-24 Single cell classification method, gene screening method and device thereof

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201110245356.8 2011-08-25
CN201110245356.8A CN102952854B (zh) 2011-08-25 2011-08-25 单细胞分类和筛选方法及其装置

Publications (1)

Publication Number Publication Date
WO2013026411A1 true WO2013026411A1 (zh) 2013-02-28

Family

ID=47745956

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2012/080578 WO2013026411A1 (zh) 2011-08-25 2012-08-24 单细胞分类方法、基因筛选方法及其装置

Country Status (5)

Country Link
US (1) US20140206006A1 (zh)
EP (1) EP2749655B2 (zh)
CN (1) CN102952854B (zh)
HK (1) HK1178210A1 (zh)
WO (1) WO2013026411A1 (zh)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111986731A (zh) * 2020-08-05 2020-11-24 广西大学 一种提高gwas致因突变定位效率的方法
WO2021042235A1 (zh) * 2019-09-02 2021-03-11 北京哲源科技有限责任公司 疾病类型自动确定方法及电子设备
CN113160886A (zh) * 2021-04-02 2021-07-23 山东大学 基于单细胞Hi-C数据的细胞类型预测系统

Families Citing this family (47)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104462869B (zh) * 2014-11-28 2017-12-26 天津诺禾致源生物信息科技有限公司 检测体细胞单核苷酸突变的方法和装置
GB2551091A (en) * 2015-04-01 2017-12-06 Toshiba Kk Genotype determination device and method
US10395759B2 (en) 2015-05-18 2019-08-27 Regeneron Pharmaceuticals, Inc. Methods and systems for copy number variant detection
CN104978497B (zh) * 2015-06-24 2019-01-29 中国科学院昆明动物研究所 一种大样本高通量生物数据关联分析方法
CN104946765A (zh) * 2015-06-25 2015-09-30 华中农业大学 基于基因组测序的体细胞突变位点挖掘方法
CN106326689A (zh) * 2015-06-25 2017-01-11 深圳华大基因科技服务有限公司 确定群体中受到选择作用的位点的方法和装置
TWI582631B (zh) * 2015-11-20 2017-05-11 財團法人資訊工業策進會 用以分析細菌菌種之定序資料的系統及其方法
CN105420374B (zh) * 2015-12-22 2019-08-16 武汉菲沙基因信息有限公司 一种诱导全能干细胞应用前期突变检测方法
CN106909806B (zh) * 2015-12-22 2019-04-09 广州华大基因医学检验所有限公司 定点检测变异的方法和装置
CN115273970A (zh) 2016-02-12 2022-11-01 瑞泽恩制药公司 用于检测异常核型的方法和系统
CN106021994B (zh) * 2016-05-13 2019-03-26 万康源(天津)基因科技有限公司 一种肿瘤突变位点筛选及互斥基因挖掘的方法
CN106021987B (zh) * 2016-05-24 2019-04-09 人和未来生物科技（长沙）有限公司 超低频突变分子标签聚类分群算法
WO2018086045A1 (zh) * 2016-11-10 2018-05-17 深圳华大基因研究院 一种对特定群中的亚群进行定量分析的方法
CN106701995B (zh) * 2017-02-20 2019-11-26 元码基因科技（北京）股份有限公司 通过单细胞转录组测序进行细胞质量控制的方法
CN106995845B (zh) * 2017-04-01 2020-05-05 中国科学院遗传与发育生物学研究所 利用三代测序平台（PacBio RS II）进行多倍体中基因等位变异挖掘的方法
CN110462063B (zh) * 2017-05-23 2023-06-23 深圳华大生命科学研究院 一种基于测序数据的变异检测方法、装置和存储介质
CN110211631B (zh) * 2018-02-07 2024-02-09 深圳先进技术研究院 一种全基因组关联分析方法、系统及电子设备
JP2019195304A (ja) 2018-05-10 2019-11-14 学校法人順天堂 画像解析方法、装置、コンピュータプログラム、及び深層学習アルゴリズムの生成方法
CN108681661B (zh) * 2018-05-16 2022-02-18 杭州迈迪科生物科技有限公司 一种伴随实验生成样本id的方法
CN108920893B (zh) * 2018-09-06 2019-04-16 南京医科大学 一种基于人工智能的颅颌面骨骼和软组织形态预测方法
CN109658981B (zh) * 2018-12-10 2022-10-04 海南大学 一种单细胞测序的数据分类方法
WO2020154885A1 (zh) * 2019-01-29 2020-08-06 北京大学 单细胞类型检测方法、装置、设备和存储介质
CN109949868B (zh) * 2019-03-01 2020-10-16 深圳乐土生物科技有限公司 基于耐受性分析的基因等级排序方法和装置
JP7381003B2 (ja) * 2019-04-26 2023-11-15 学校法人順天堂 疾患解析を支援する方法、装置、及びコンピュータプログラム、並びにコンピュータアルゴリズムを訓練する方法、装置、及びプログラム
CN110504005A (zh) * 2019-08-27 2019-11-26 上海其明信息技术有限公司 数据处理方法
CN110797087B (zh) * 2019-10-17 2020-11-03 南京医基云医疗数据研究院有限公司 测序序列处理方法及装置、存储介质、电子设备
CN110827921B (zh) * 2019-11-12 2022-06-14 玉林师范学院 一种单细胞聚类方法、装置、电子设备及存储介质
CN111105844B (zh) * 2019-11-22 2023-06-06 广州金域医学检验集团股份有限公司 体细胞变异分类方法、装置、设备及可读存储介质
US20230028790A1 (en) * 2019-11-29 2023-01-26 GC Genome Corporation Artificial intelligence-based chromosomal abnormality detection method
CN110993027B (zh) * 2019-12-18 2022-10-11 武汉大学 一种高效克隆植物性状相关突变基因的方法
CN111091867B (zh) * 2019-12-18 2021-11-09 中国科学院大学 基因变异位点筛选方法及系统
CN111681710B (zh) * 2020-06-03 2021-08-27 中国人民解放军军事科学院军事医学研究院 基于基因表达特征的细胞分类方法、装置和电子设备
CN111899790A (zh) * 2020-08-17 2020-11-06 天津诺禾医学检验所有限公司 测序数据的处理方法及装置
CN112289379B (zh) * 2020-10-15 2022-11-22 天津诺禾致源生物信息科技有限公司 细胞类型的确定方法、装置、存储介质及电子装置
CN112270957B (zh) * 2020-10-19 2023-11-07 西安邮电大学 高阶snp致病组合数据检测方法、系统、计算机设备
CN112837754B (zh) * 2020-12-25 2022-10-28 北京百奥智汇科技有限公司 一种基于特征基因的单细胞自动分类方法和装置
CN113223611B (zh) * 2021-04-16 2024-02-20 杭州电子科技大学 一种基于全基因组序列snp数据的快速特征选择方法
CN113674800B (zh) * 2021-08-25 2022-02-08 中国农业科学院蔬菜花卉研究所 基于单细胞转录组测序数据的细胞聚类方法
WO2023157933A1 (ja) * 2022-02-18 2023-08-24 国立研究開発法人理化学研究所 細胞の系統解析を行う方法
CN115394358B (zh) * 2022-08-31 2023-05-12 西安理工大学 基于深度学习的单细胞测序基因表达数据插补方法和系统
GB2622371A (en) * 2022-09-13 2024-03-20 Agecurve Ltd Cell tree rings: Method and cell lineage tree based aging timer for calculating biological age of biological sample
CN115424668B (zh) * 2022-11-02 2023-03-24 杭州联川基因诊断技术有限公司 一种单细胞转录组数据可用性分析方法、介质及设备
CN116230082B (zh) * 2022-12-06 2024-05-14 序科码医学检验实验室(广州)有限公司 基于样品基因型的数据拆分的无标记多样品混合单细胞测序方法
CN115966259B (zh) * 2022-12-26 2023-10-13 南京普恩瑞生物科技有限公司 一种基于逻辑回归建模的样本同源性检测校验方法及系统
CN116486913B (zh) * 2023-05-23 2023-10-03 浙江大学 基于单细胞测序从头预测调控突变的系统、设备和介质
CN117423382B (zh) * 2023-10-21 2024-05-10 云准医药科技(广州)有限公司 一种基于SNP多态性的单细胞barcode身份识别方法
CN117637020B (zh) * 2024-01-25 2024-04-30 鲁东大学 一种基于深度学习的四倍体牡蛎全基因组snp分型方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101914628A (zh) * 2010-09-02 2010-12-15 深圳华大基因科技有限公司 检测基因组目标区域多态性位点的方法及 系统
CN101956006A (zh) * 2010-08-27 2011-01-26 公安部物证鉴定中心 获得人种特异性位点的方法和人种推断系统及其应用
CN102061526A (zh) * 2010-11-23 2011-05-18 深圳华大基因科技有限公司 一种DNA文库及其制备方法、以及一种检测SNPs的方法和装置

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101956006A (zh) * 2010-08-27 2011-01-26 公安部物证鉴定中心 获得人种特异性位点的方法和人种推断系统及其应用
CN101914628A (zh) * 2010-09-02 2010-12-15 深圳华大基因科技有限公司 检测基因组目标区域多态性位点的方法及 系统
CN102061526A (zh) * 2010-11-23 2011-05-18 深圳华大基因科技有限公司 一种DNA文库及其制备方法、以及一种检测SNPs的方法和装置

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
"Sequencing of 50 Human Exomes Reveals Adaptation to High Altitude", SCIENCE, vol. 329, 2 July 2010 (2010-07-02), pages 75 - 78
J. SAMBROOK: "Molecular Cloning: A Laboratory Manual, 3rd edition,", SCIENCE PRESS
See also references of EP2749655A4

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
WO2021042235A1 (zh) * 2019-09-02 2021-03-11 北京哲源科技有限责任公司 疾病类型自动确定方法及电子设备
CN112930573A (zh) * 2019-09-02 2021-06-08 北京哲源科技有限责任公司 疾病类型自动确定方法及电子设备
CN111986731A (zh) * 2020-08-05 2020-11-24 广西大学 一种提高gwas致因突变定位效率的方法
CN111986731B (zh) * 2020-08-05 2023-08-11 广西大学 一种提高gwas致因突变定位效率的方法
CN113160886A (zh) * 2021-04-02 2021-07-23 山东大学 基于单细胞Hi-C数据的细胞类型预测系统

Also Published As

Publication number Publication date
EP2749655B2 (en) 2022-04-13
CN102952854B (zh) 2015-01-14
EP2749655A1 (en) 2014-07-02
US20140206006A1 (en) 2014-07-24
HK1178210A1 (zh) 2013-09-06
EP2749655B1 (en) 2019-07-03
EP2749655A4 (en) 2015-01-28
CN102952854A (zh) 2013-03-06

Similar Documents

Publication Publication Date Title
WO2013026411A1 (zh) 单细胞分类方法、基因筛选方法及其装置
US10947595B2 (en) Nucleic acids and methods for detecting chromosomal abnormalities
US20200251180A1 (en) Resolving genome fractions using polymorphism counts
TWI661049B (zh) 使用不含細胞之ｄｎａ片段大小以測定複製數變異之方法
CN110800063B (zh) 使用无细胞dna片段大小检测肿瘤相关变体
Fan et al. Non-invasive prenatal measurement of the fetal genome
Clark et al. Performance comparison of exome DNA sequencing technologies
CN113366122B (zh) 游离dna末端特征
KR20210003094A (ko) 잔류 질환의 검출을 위한 시스템 및 방법
CN107771221A (zh) 用于癌症筛查和胎儿分析的突变检测
CN104884633B (zh) 通过测序少量遗传物质的高通量基因分型
US20210090687A1 (en) Methods of quality control using single-nucleotide polymorphisms in pre-implantation genetic screening
JP2022537445A (ja) 精子提供者、卵母細胞提供者、及びそれぞれの受胎産物の間の遺伝的関係を決定するためのシステム、コンピュータプログラム製品及び方法
IT202100024101A1 (it) Metodo per analizzare il grado di similarita' di almeno due campioni utilizzando amplificazione deterministica dell'intero genoma mediante siti di restrizione (drs-wga)
WO2024192076A1 (en) Sample barcode in multiplex sample sequencing
Shani et al. Cell-free DNA: Comparison of Technologies
CN118460721A (zh) 一种原发性中枢神经系统淋巴瘤分子分型试剂盒及其应用

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 12825673

Country of ref document: EP

Kind code of ref document: A1

WWE Wipo information: entry into national phase

Ref document number: 14239650

Country of ref document: US

NENP Non-entry into the national phase

Ref country code: DE