Method for detecting genetic variation in highly homologous sequences by independent alignment and pairing of sequence reads
By performing peer-sequencing and top-pairing alignments on the genome, the problem of inaccurate detection of highly homologous gene regional variants in the prior art is solved, and rapid and accurate gene variant detection is achieved.
Patent Information
- Application Number
- JP2025019414
- Authority / Receiving Office
- JP · JP
- Patent Type
- Applications
- Current Assignee / Owner
- Priority Date
- 2018-09-12
- Filing Date
- 2025-02-07
- Publication Date
- 2025-05-15
AI Technical Summary
The prior art is difficult to accurately detect genetic mutations, especially in highly homologous gene regions, resulting in inaccurate detection results.
By performing peer-sequencing of the genome, the first and second read sequences are obtained and aligned to the reference genome respectively to generate top paired alignments, and then gene mutations are detected.
This method can quickly and efficiently detect genetic mutations, especially in highly homologous regions, improving the accuracy and efficiency of detection.
Smart Images

Figure 2025076479000001_ABST
Abstract
Description
[Technical field]
[0001] CROSS-REFERENCE TO RELATED APPLICATIONS
[0001] This application claims priority to U.S. Provisional Application No. 62 / 711,454, filed July 27, 2018, and U.S. Provisional Application No. 62 / 730,479, filed September 12, 2018, each of which is incorporated herein in its entirety, including all tables, figures, and claims.
[0002]
[0002] The following disclosure relates generally to determining genetic variation, and more particularly to determining genetic variation in highly homologous regions of interest in a genome, e.g., in genomic regions that contain genes and pseudogenes. [Background technology]
[0003]
[0003] Individual genomic variants inherited through the germline account for approximately 5% to 10% of cancers [1-3]. This hereditary component can increase the risk of malignancies across a range of tissues [4, 5] (e.g., breast, colorectal, pancreatic, and prostate) and has been associated with pathogenic variants in over 100 genes [6]. To assess a patient's risk for such cancers, hereditary cancer screening (HSC) typically uses targeted next-generation sequencing (NGS) to detect associated variants in coding regions and select non-coding regions in multigene testing panels.
[0004]
[0004] In most genomic regions surveyed by HSC panels, NGS alone is sufficient to obtain high sensitivity and specificity [7, 8], and high accuracy is important for HSC because the results of the test prompt patients to change their clinical management decisions [9, 10]. However, in a few regions, standard NGS strategies that use hybridization to capture and sequence short DNA fragments can only identify genotypes inaccurately. Genes with particular challenges often have sequences homologous to other parts of the genome (e.g., pseudogenes) that are captured and sequenced together with the gene itself, complicating alignment and identification of gene-specific variants.
[0005]
[0005] Thus, there remains a need for improved methods for detecting genetic variations in homologous regions of the genome. Summary of the Invention
[0006]
[0006] Current techniques that allow for the determination of genotypes for highly homologous genes and corresponding homologs are time-consuming, laborious, and expensive, making them unsuitable for widespread clinical use.
[0007]
[0007] The disclosed method can be implemented in an affordable and high-throughput manner, thus saving considerable time, effort and cost. Moreover, the method overcomes the challenge of elucidating the structure / copy number / genotype in regions where the unique alignment of NGS reads to genes or their homologs is compromised.
[0008]
[0008] In one aspect, provided herein is a method for determining an individual's genomic structure (i.e., genotype) with respect to a gene of interest, wherein the gene of interest has a highly homologous homolog, e.g., a pseudogene.
[0009]
[0009] In one embodiment, a method for detecting genetic variations in a genome of a subject, the genome comprising a first and a second region of interest that are highly homologous, the method comprising: (a) obtaining sequence reads by paired-end sequencing from a number of sites of interest in the first and second regions of interest, the sequence reads comprising a first read and a second read obtained at each site of interest; (b) aligning the sequence reads to a reference genome, the first read and the second read being separately aligned to the reference genome, and an aligner generating a number of possible alignments for each of the first read and the second read; (c) identifying the first read and the second read that align to the first region of interest; and (d) pairing the first read and the second read from the reads identified in step (c), thereby generating a top paired alignment. A method is provided, comprising: (a) generating a top pair alignment of the first and second reads; and (b) detecting genetic variations in the top pair alignment generated in step (d). In another embodiment, the method includes, prior to step (b), aligning the first and second reads to a reference genome, where the aligner generates a best possible paired end alignment for each pair of the first and second reads to the first or second region of the object, and only the paired end reads associated with the top alignment score for the first or second region of the object are aligned separately in step (b). In one embodiment, the reference genome does not include masked or modified portions of the first or second homologous region of the object. In one embodiment, the method is implemented by a computer.
[0010]
[0010] In one embodiment, a method is provided for detecting genetic variations in a subject's genome, the genome comprising highly homologous first and second regions of interest, the method comprising the step of obtaining sequence reads by paired-end sequencing from multiple sites of interest in the first and second regions of interest, the sequence reads comprising a first read and a second read obtained at each site of interest, the sequence reads being obtained by direct targeted sequencing (DST) of the multiple sites of interest, and the first read comprising a genomic sequence and the second read comprising a probe sequence read associated with the sites of interest.
[0011]
[0011] In one embodiment, a method is provided for detecting genetic variations in a subject's genome, the genome comprising a first and second region of interest that are highly homologous to the target, the method comprising: (a) obtaining sequence reads by paired-end sequencing from a number of sites of interest in the first and second regions of the target, the sequence reads comprising a first read and a second read obtained at each site of the target; (b) aligning the sequence reads to a reference genome, the first read and the second read being separately aligned to the reference genome, and an aligner generating a number of possible alignments for each of the first read and the second read; (c) identifying the first read and the second read that align to the first region of the target; (d) pairing the first read and the second read from the reads identified in step (c), thereby generating a top pair alignment; and (e) detecting genetic variations in the top pair alignment generated in step (d). In one embodiment, the sequence reads are aligned using the Burrows-Wheeler Aligner (BWA) algorithm. In one embodiment, the aligner only issues alignments that meet a minimum alignment score for the first and second regions of the object. In one embodiment, the first and second reads are paired to produce a top pair alignment only if the alignments of the first and second reads to the first region of the object are within a certain number of bases of each other. In one embodiment, the first and second reads are paired to produce a top pair alignment only if the alignments of the first and second reads to the first region of the object are within a certain number of bases of each other. The reads are paired to produce a top pair alignment only if the alignment of the first read and the second read to the first region of interest is within about 100 bp, about 200 bp, about 200 bp, about 300 bp, about 400 bp, about 500 bp, about 600 bp, about 700 bp, about 800 bp, about 900 bp, about 1000 bp, about 1100 bp, about 1200 bp, about 1300 bp, about 1400 bp, about 1500 bp, or more than 1500 bp.
[0012]
[0012] In one embodiment, a method is provided for detecting genetic variations in a genome of a subject, the genome comprising a first and second region of interest that are highly homologous to the target, the method comprising: (a) obtaining sequence reads by paired-end sequencing from a number of sites of interest in the first and second regions of the target, the sequence reads comprising a first read and a second read obtained at each site of the target; (b) aligning the sequence reads to a reference genome, the first read and the second read being separately aligned to the reference genome, and an aligner generating a number of possible alignments for each of the first read and the second read; (c) identifying the first read and the second read that align to the first region of the target; (d) pairing the first read and the second read from the reads identified in step (c), thereby generating a top pair alignment; and (e) detecting genetic variations in the top pair alignment generated in step (d). In one embodiment, the method includes, in step (d), generating a number of pairwise alignments, calculating an alignment score for each of the number of pairwise alignments, and identifying a top pairwise alignment having the highest alignment score. In one embodiment, the top pairwise alignment in step (d) is selected as the one having the smallest template length.
[0013]
[0013] In one embodiment, a method is provided for detecting genetic variations in a subject's genome, the genome comprising a first and second region of interest that are highly homologous to the target, the method comprising: (a) obtaining sequence reads by paired-end sequencing from a number of sites of interest in the first and second regions of interest, the sequence reads comprising a first read and a second read obtained at each site of the target; (b) aligning the sequence reads to a reference genome, the first read and the second read being separately aligned to the reference genome, and an aligner generating a number of possible alignments for each of the first read and the second read; (c) identifying the first read and the second read that align to the first region of the target; (d) pairing the first read and the second read from the reads identified in step (c), thereby generating a top pair alignment; and (e) detecting genetic variations in the top pair alignment generated in step (d). In one embodiment, the genetic variation comprises a SNP, an indel, an inversion, and / or a CNV. In one embodiment, the detecting step in step (e) comprises calling the SNP, the indel, the inversion, and / or the CNV. In one embodiment, the detecting step in step (e) comprises using a Hidden Markov Model (HMM) caller to determine copy number. In one embodiment, the detecting step in step (e) is based on an expected ploidy of 2. In one embodiment, the detecting step in step (e) is based on an expected ploidy of 4. In one embodiment, when a genetic variation is detected in step (e), a portion of the genome of the subject is amplified by long-range PCR and assayed by multiplex ligation-dependent probe amplification (MLPA). In one embodiment, when a genetic variation is detected in step (e), a portion of the first region of interest is amplified by long-range PCR and the product or a portion thereof is subjected to Sanger sequencing or NGS. In one embodiment, if a genetic mutation is detected in step (e), the subject's genomic DNA is assayed by multiplex ligation-dependent probe amplification (MLPA).
[0014]
[0014] In one embodiment, a method is provided for detecting genetic variations in a subject's genome, the genome comprising a first and second region of interest that are highly homologous to the target, the method comprising: (a) obtaining sequence reads by paired-end sequencing from a number of sites of interest in the first and second regions of the target, the sequence reads comprising a first read and a second read obtained at each site of the target; (b) aligning the sequence reads to a reference genome, the first read and the second read being separately aligned to the reference genome, and an aligner generating a number of possible alignments for each of the first read and the second read; (c) identifying the first read and the second read that align to the first region of the target; (d) pairing the first read and the second read from the reads identified in step (c), thereby generating a top pair alignment; and (e) detecting genetic variations in the top pair alignment generated in step (d). In one embodiment, the sequence reads are 30-50 bp or 100-200 bp in length. In one embodiment, the highly homologous first and second regions of interest are at least 80%, at least 81%, at least 82%, at least 83%, at least 84%, at least 85%, at least 86%, at least 87%, at least 88%, at least 89%, at least 90%, at least 91%, at least 92%, at least 93%, at least 94%, at least 95%, at least 96%, at least 97%, at least 98%, at least 99%, or more than 99% identical. In one embodiment, the sequence reads are obtained from one or more exons in the first and / or second region of interest. In one embodiment, the sequence reads are obtained from one or more introns in the first and / or second region of interest. In one embodiment, sequence reads are obtained from one or more exons and introns within the first and / or second region of interest.In one embodiment, the sequence reads are obtained from one or more exons and introns in the first and / or second region of interest, where the introns are located near the exons. In one embodiment, the sequence reads are obtained from one or more clinically addressable regions associated with the first and / or second region of interest. In one embodiment, the first region of interest comprises a gene and the second region of interest comprises a pseudogene. In one embodiment, the first region of interest comprises a pseudogene and the second region of interest comprises a gene. In one embodiment, the first region of interest comprises two alleles. In one embodiment, the second region of interest comprises two alleles. In one embodiment, the gene is PMS2. In one embodiment, the pseudogene is PMS2CL. In one embodiment, the multiple sites of interest are located within an exon of PMS2 and an exon of another portion of the genome of the subject. In one embodiment, the multiple sites of interest are located within an exon of PMS2 and an exon of PMS2CL. In one embodiment, multiple sites of interest are present within exons 11, 12, 13, 14, and / or 15 of PMS2 and exons 2, 3, 4, 5, and / or 6 of PMS2CL. In one embodiment, the subject is a human and the sequence reads are aligned to a human reference genome.
[0015]
[0015] In one embodiment, a method for detecting genetic variations in a genome of a subject, the genome comprising first and second highly homologous regions of interest, the method comprising: (a) obtaining sequence reads from multiple sites of interest in the first and second regions of interest by paired-end sequencing, the sequence reads comprising a first read and a second read obtained at each site of the interest; and (b) obtaining a reference A method is provided that includes aligning sequence reads to a genome, where a first read and a second read are aligned to a reference genome separately, and an aligner generates multiple possible alignments for each of the first read and the second read; (c) identifying the first read and the second read that align to a first region of interest; (d) pairing the first read and the second read from the reads identified in step (c), thereby generating a top pair alignment; and (e) detecting genetic variations in the top pair alignment generated in step (d). In one embodiment, a non-transitory computer-readable storage medium is provided that includes computer-executable instructions for performing the method described herein. In one embodiment, a system is provided that includes (a) one or more processors, (b) a memory, and (c) one or more programs, where the one or more programs are stored in the memory and configured to be executed by the one or more processors, and the one or more programs include instructions for performing the method described herein.
[0016] In one embodiment, a computer system is provided that is configured to execute instructions for performing the methods described herein.
[0017]
[0017] Other objects, features and advantages of the present invention will become apparent from the following detailed description. However, it should be understood that the detailed description and specific examples, while indicating preferred embodiments of the present invention, are given by way of illustration only, since various changes and modifications within the scope and spirit of the present invention will become apparent to those skilled in the art from this detailed description. [Brief description of the drawings]
[0018] [Figure 1A]
[0018] Figures 1A-1D show the LR-PCR strategy for constructing a dataset of natural genetic variations in PMS2 and PMS2CL. Figure 1A: Short reads from NGS hybrid-capture data originating from genes (blue) and pseudogenes (red) align to both genes and pseudogenes due to high homology. [Figure 1B] Figures 1A-1D show the LR-PCR strategy for constructing a dataset of natural genetic variation in PMS2 and PMS2CL. Figure 1B: Using LR-PCR that is specific for genes or pseudogenes, followed by fragmentation and barcoding (Figure 1B), the resulting NGS short reads can be assigned to genes or pseudogenes (Figure 1C). [Figure 1C] Figures 1A-1D show the LR-PCR strategy for constructing a dataset of natural genetic variation in PMS2 and PMS2CL. Figure 1C: Using LR-PCR that is specific for genes or pseudogenes, followed by fragmentation and barcoding (Figure 1B), the resulting NGS short reads can be assigned to genes or pseudogenes (Figure 1C). [Figure 1D] Figures 1A-1D show the LR-PCR strategy for constructing a dataset of natural genetic variations in PMS2 and PMS2CL. Figure 1D: Percent identity between the gene and pseudogene for PMS2 exons 11-15 based on the hg19 reference genome (grey) and after taking into account the natural genetic variations from the LR-PCR samples (black). [Figure 2A]
[0019] Figures 2A-2B show the reflex workflow for variant identification in the final exons of PMS2. Figure 2A: Overview of the sequencing and analysis workflow for the five final exons of PMS2. Colored nodes correspond to the boxes in Figure 2B. [Figure 2B]Figures 2A-2B show the reflex workflow for variant identification in the final exon of PMS2. Figure 2B: Details corresponding to steps of the workflow in Figure 2A; details for each box are described in Methods and Results. "Not reported" means that the variant does not appear in the patient's report. "Reflex" means that the sample is sent for LR-PCR-based disambiguation to determine whether the variant localizes to a gene or pseudogene. [Figure 3A]
[0020] Figures 3A-3C show that hybrid-capture and LR-PCR are responsive to SNVs and indels. Figure 3A: Hypothetical example illustrating a correspondence table for comparison of hybrid capture and LR-PCR data. All examples assume that the reference base is A and the alternative ("alt") base is T. (i) Example of a true positive (dark blue) where the alt allele is present in PMS2CL. (ii) Example of allowed dosage error (light blue) where PMS2CL is homozygous for the alt allele, but hybrid capture calls only one alt allele instead of two. (iii) Example of a false positive (light orange) where only hybrid capture detected the alt allele. (iv) Example of a false negative (dark orange) where the alt allele in PMS2CL was not captured by hybrid capture. The shaded matrix on the right shows cells representing true positives, allowed dosage error, false positives, and false negatives. Numbers on the axis indicate the total number of alt alleles in either the hybrid capture data or the PMS2 / PMS2CL LR-PCR data. [Figure 3B] Figures 3A-3C show hybrid-capture and LR-PCR correspond to SNVs and indels. Figure 3B: Diploid SNVs and indels correspond to exon 11 of PMS2. Numbers on axis indicate number of alt alleles with 0 equal to 0 / 0, 1 equal to 0 / 1, and 2 equal to 1 / 1. 95% confidence intervals in brackets. [Figure 3C]Figures 3A-3C show that hybrid-capture and LR-PCR correspond to SNVs and indels. Figure 3C: Four copies of SNVs and indels correspond to exons 12-15 of PMS2 / PMS2CL as described in Figure 3A. [Figure 4A]
[0021] Figures 4A-4B show that simulated indels increase confidence in indel sensitivity. Figure 4A: Schematic of simulating indels in tetraploids by combining sequencing data from two diploid samples. [Figure 4B] Figures 4A-4B show that simulated indels increase confidence in indel sensitivity. Figure 4B: Simulation results for tetraploid indels in the same format as Figure 3A. [Figure 5A]
[0022] Figures 5A-5D show hybrid capture, LR-PCR, and MLPA correspond to CNVs. Figure 5A: All CNVs called in hybrid capture data and corresponding orthogonal validation data. [Figure 5B] Figures 5A-5D show hybrid capture, LR-PCR, and MLPA corresponding to CNVs. Figure 5B: Hybrid capture data for a patient sample with an exon 13-14 deletion shows copy number estimates across loci (bins). Grey areas indicate the four final exons of PMS2. White areas indicate introns. Yellow boxes indicate regions of CNV calls. [Figure 5C] Figures 5A-5D show hybrid capture, LR-PCR, and MLPA corresponding to the CNV. Figure 5C: MLPA data for a patient sample with a deletion of exons 13-14. The PMS2-specific MLPA probe (solid blue), PMS2CL-specific MLPA probe (solid red), and PMS2 / PMS2CL degenerate MLPA probe (blue and red stripes) show a deletion in exons 13-14 of PMS2CL. [Figure 5D]Figures 5A-5D show hybrid capture, LR-PCR, and MLPA corresponding to CNVs. Figure 5D: LR-PCR data for exon 13-14 deletion sample showing copy number estimates across loci (bins) for PMS2 (blue, top) and PMS2CL (red, bottom). Grey areas indicate exons 11-15 of PMS2, white areas indicate introns as in Figure 5B. [Figure 6]
[0023] Figure 6 shows the orthogonal data sets used to construct the hybrid capture assay. As shown, Figure 6 is a diagram that demonstrates the assays, data sets, algorithms, and analyses used to construct the hybrid capture assay for the five final exons of PMS2. The Coriell sample (1b) can be used by other researchers without repeating the LR-PCR provided in accession number PRJEB27948. Genomic DNA (gDNA). [Figure 7]
[0024] Figures 7A-7C show that the PMS2 exon 11-15 reference genotypes (from Polaris and GIAB) are not concordant with the PMS2 LR-PCR. Figure 7A: Concordance between LR-PCR variant calls and Polaris variant cells. Figure 7B: Concordance between LR-PCR variant cells and the GIAB multi-sample call set (including high confidence and filtered variant cells) for all five GIAB samples. Figure 7C: Concordance between LR-PCR variant calls and the available 10x Genomics haplotype call set for four GIAB samples. [Figure 8A]
[0025] Figures 8A-8B show that the RNA data corroborates the hybrid capture and LR-PCR data. Figure 8A: Concordance between hybrid capture data and RT-PCR for PMS2 and PMS2CL. [Figure 8B] Figures 8A-8B show that the RNA data corroborates the hybrid capture and LR-PCR data. Figure 8B: Concordance between hybrid capture data and LR-PCR for PMS2 and PMS2CL. [Figure 9]
[0026] FIG. 9 is a chart illustrating an embodiment of a method described herein that involves an "ambiguous alignment" of a first DTS read and a second DTS read from a region of interest. [Figure 10]
[0027] FIG. 10 illustrates an example system and environment in which various embodiments of the invention may operate. [Figure 11]
[0028] FIG. 11 is a diagram illustrating an exemplary computing system. DETAILED DESCRIPTION OF THE PREFERRED EMBODIMENTS
[0019]
[0029] The file of this patent contains at least one drawing in color. Copies of this patent or patent publication with color drawing(s) will be provided by the Office upon application and payment of the necessary fee.
[0020]
[0030] The present invention will now be described in detail, by way of reference only, by using the following definitions and examples. All patents and publications referred to herein are expressly incorporated by reference, including all sequences disclosed within such patents and publications.
[0021]
[0031] Unless otherwise defined herein, all technical and scientific terms used herein have the same meaning as commonly understood by those skilled in the art to which this invention belongs. Singleton et al., Dictionary of Microbiology and Molecular Biology, 2nd ed., John Wiley and Sons, New York (1994), and Hale and Marham, The Harper Collins Dictionary of Biology, Harper Perennial, NY (1991) provide those skilled in the art with a general dictionary for many of the terms used in this invention. Although any methods and materials similar or equivalent to those described herein can be used in the practice or testing of this invention, preferred methods and materials are described. In particular, experts are directed to Sambrook et al., 1989, and Ausubel FM et al., 1993 for definitions and terms of the art. It should be understood that this invention is not limited to the specific methodology, protocols, and reagents described, as they may vary.
[0022]
[0032] Numeric ranges are inclusive of the numbers that define the range. The term "about" is used herein to mean plus or minus ten percent (10%) of a value. For example, "about 100" refers to any number between 90 and 110.
[0023]
[0033] Unless otherwise indicated, nucleic acids are written left to right in 5' to 3' orientation; amino acid sequences are written left to right in amino to carboxy orientation, respectively.
[0024]
[0034] The headings provided herein are not limitations of the various aspects or embodiments of the invention that may have by reference to the specification as a whole. Accordingly, the terms defined immediately below are more fully defined by reference to the specification as a whole.
[0035] Supplemental data, including any referenced tables (e.g., Table S1, Table S2, etc.), will be available upon request. A version of the scientific article relating to this patent application is provided as an appendix hereto.
[0025] I. Definition
[0036] As used herein, "purified" and its derivatives mean that a molecule is present in a sample at a concentration of at least 90%, 95%, or at least 98% by weight of the sample in which it is contained.
[0026]
[0037] The term "isolated" and its derivatives, as used herein, refers to a molecule that is separated from at least one other molecule with which it is normally associated, e.g., in a natural environment. An isolated nucleic acid molecule includes a nucleic acid molecule that is originally contained within a cell that normally expresses the nucleic acid molecule, but where the nucleic acid molecule is present extrachromosomally or at a chromosomal location that is different from its natural chromosomal location.
[0027]
[0038] The term "% identity" and its derivatives are used interchangeably herein with the term "% homology" and its derivatives to refer to the level of identity of a nucleic acid or amino acid sequence between another nucleic acid sequence or any other polypeptide, or the amino acid sequence of a polypeptide, to which the sequence is aligned using a sequence alignment program, e.g., using the Basic Local Alignment Search Tool algorithm. In the case of nucleic acids, the term also applies to intronic and / or intergenic regions.
[0028]
[0039] For example, as used herein, 80% homology means the same as 80% sequence identity as determined by a defined algorithm, and thus a homologue or highly homologous sequence of a given sequence has a higher percentage sequence identity than 80% over the length of the given sequence. Exemplary levels of sequence identity include, but are not limited to, 80, 85, 90, 95, 98% or higher percentage sequence identity to a given sequence, for example, a coding sequence for any one of the polypeptides of the present invention, as described.
[0029]
[0040] As used herein, "highly homologous" and its derivatives mean that the percent homology or percent identity between at least two different nucleotide sequences is greater than 70%. Sequences are referred to as "highly homologous" if their sequence identity is greater than 70% over comparable lengths.
[0030]
[0041] Exemplary computer programs that can be used to determine identity between two sequences include, but are not limited to, the BLAST suite of programs, such as BLASTN, BLASTX, and TBLASTX, BLASTP, and TB Examples of such protocols include BLAST, LASTN, and BLAST, which are publicly available on the internet. See also Altschul et al., 1990 and Altschul et al., 1997.
[0031]
[0042] Sequence searches are typically performed using the BLASTN program, where a given nucleic acid sequence is evaluated against the DNA sequences in GenBank and other public databases. The BLASTX program is preferred for searching nucleic acid sequences translated in all reading frames against the amino acid sequences in GenBank and other public databases. Both BLASTN and BLASTX are performed using default parameters, with an open gap penalty of 11.0 and an extended gap penalty of 1.0, and utilize the BLOSUM-62 matrix. (See, e.g., Altschul, SF et al., Nucleic Acids Res. 25:3389-3402, 1997).
[0032]
[0043] Preferred alignment of selected sequences to determine the "% identity" between two or more sequences is performed, for example, using the CLUSTAL-W program in MacVector version 13.0.7, operated with default parameters, including an open gap penalty of 10.0, an extended gap penalty of 0.1, and a BLOSUM 30 similarity matrix.
[0033]
[0044] A "sequence read" and its derivatives range from 30 nt to 400 nt, 50 nt to 250 nt, 50 nt to 150 nt, or 100 nt to 200 nt in a nucleotide sequence.
[0034]
[0045] The term "mutation," as used herein, refers to both natural and inherited sequence variations, including, but not limited to, variations between individuals, or between an individual's sequence and a reference sequence. Exemplary mutations include, but are not limited to, SNPs, indels (insertion or deletion variants), copy number variants, inversions, translocations, chromosomal fusions, and the like.
[0035]
[0046] The term "small nucleotide polymorphism" or "SNP" and its derivatives refer to single nucleotide variants (SNVs), multi-nucleotide variants (MNVs), or indel variants of about 100 base pairs or less.
[0036]
[0047] The term "homolog" and its derivatives, as used herein, refers to a nucleotide sequence that is identical or nearly identical to a nucleotide sequence located elsewhere in the genome of interest. A homolog is highly homologous to a nucleotide sequence located elsewhere in the genome of interest. A homolog can be either a "pseudogene" which is another gene, or a segment of sequence that is not part of a gene.
[0037]
[0048] "Pseudogene" and its derivatives, as used herein, is a DNA sequence that closely resembles a gene in the DNA sequence, but has at least one change that renders the gene dysfunctional. The change may be a single residue mutation. The change may result in a splice variant. The change may result in premature termination of translation. A pseudogene is dysfunctional relative to a functional gene. A pseudogene is characterized by a combination of homology to a known gene (i.e., a gene of interest) and non-functionality.
[0038]
[0049] The number of pseudogenes for a gene is not limited to those enumerated herein. Pseudogenes are increasingly being recognized. Therefore, those skilled in the art can easily identify pseudogenes based on sequence homology or by using, for example, GeneCards (genecards.org), pseudogenes, etc. Curated databases such as dogenes.org can be consulted to determine whether a sequence is a pseudogene.
[0039]
[0050] As used herein, "gene of interest" and its derivatives are genes for which it is desired to determine the genotype. Overall, a gene of interest has two functional copies, with two chromosomes each carrying a copy of the gene of interest. The terms "gene of interest" and "gene" can be used interchangeably herein.
[0040]
[0051] As used herein, a "region of interest" and its derivatives may be any region within the genome of a subject. As used herein, a region of interest is a sequence that is highly homologous in the genome of a subject as a whole.
[0041] II. Process
[0052] The samples in which polynucleotides are analyzed by the methods described herein may be derived from multiple samples from the same individual, multiple samples from different individuals, or combinations thereof. In some embodiments, the sample contains multiple polynucleotides from a single individual. In some embodiments, the sample contains multiple polynucleotides from two or more individuals. For example, the sample is derived from a pregnant woman and contains polynucleotides from the pregnant woman and her fetus. An individual is any organism or part thereof from which polynucleotides can be derived, non-limiting examples of which include plants, animals, fungi, protists, Monera organisms, viruses, mitochondria, and chloroplasts. Sample polynucleotides may be isolated from a subject, such as a cell sample, tissue sample, fluid sample, or organ sample derived therefrom (or cell culture derived from any of these), including cultured cell lines, biopsies, blood samples, buccal swabs, fluid samples containing cells (e.g., saliva), and the like. The subject may be an animal, including, but not limited to, cows, pigs, mice, rats, chickens, cats, dogs, and the like, and is typically a mammal, such as a human. The sample may be derived artificially, such as by chemical synthesis. In some embodiments, the sample comprises DNA. In some embodiments, the sample comprises cell-free DNA extracted from the plasma of a subject. In some embodiments, the sample comprises genomic DNA. In some embodiments, the sample comprises mitochondrial DNA, chloroplast DNA, plasmid DNA, bacterial artificial chromosomes, yeast artificial chromosomes, oligonucleotide tags, polynucleotides from organisms (e.g., bacteria, viruses, or fungi) other than the subject from which the sample is obtained, or combinations thereof. In some embodiments, the extracted nucleic acid comprises cell-free DNA from maternal plasma of a pregnant woman.
[0042]
[0053] Methods for nucleic acid extraction and purification are well known in the art. For example, nucleic acid can be purified by organic extraction, including phenol, phenol / chloroform / isoamyl alcohol, or similar preparations, including TRIzol and TriReagent. Other non-limiting examples of extraction techniques include: (1) organic extraction followed by ethanol precipitation, for example, using phenol / chloroform organic reagent (Ausubel et al., 1993), with or without the use of an automated nucleic acid extractor, for example, Model 341 DNA Extractor available from Applied Biosystems (Foster City, Calif.); (2) solid phase adsorption method (U.S. Pat. No. 5,234,809; Walsh et al., 1991); and (3) salt-induced nucleic acid precipitation, such as the precipitation method typically referred to as "salting out" method (Miller et al., (1988)). Another example of nucleic acid isolation and / or purification involves the use of magnetic particles, where the nucleic acid can be specifically or non-specifically bound to the magnetic particles, followed by isolating the beads using a magnet, washing, and eluting the nucleic acid from the beads (see, e.g., U.S. Pat. No. 5,705,628). In some embodiments, the above isolation methods may be preceded by an enzymatic digestion step, e.g., digestion with proteinase K, or other similar proteases, which helps to remove unwanted proteins from the sample. See, e.g., U.S. Pat. No. 7,001,724. In a preferred embodiment, the extracted DNA comprises the genome of the subject.
[0043]
[0054] In some embodiments, a library (e.g., a DNA library) comprising a plurality of nucleic acid molecules is prepared from the extracted nucleic acid. In some embodiments, the nucleic acids in the plurality of nucleic acid molecules comprise incorporated oligonucleotides, which may include molecular barcodes and / or one or more adaptor oligonucleotides (also referred to as "adapters").
[0044]
[0055] In some embodiments, a portion of the extracted nucleic acid is amplified by a primer extension reaction using any suitable combination of primers and DNA polymerase, including, for example, but not limited to, polymerase chain reaction (PCR), reverse transcription, and combinations thereof. When the template for the primer extension reaction is RNA, the reverse transcription product is referred to as complementary DNA (cDNA). Primers useful in primer extension reactions may include one or more target-specific sequences, random sequences, partially random sequences, and combinations thereof. Suitable reaction conditions for primer extension reactions are known in the art. In some embodiments, the extracted DNA is amplified by long-range PCR (LR-PCR) using specific primers, for example, gene-specific primers.
[0045]
[0056] The extracted nucleic acid is sequenced. Methods for sequencing nucleic acids are well known in the art. In one embodiment, the extracted nucleic acid is sequenced by Sanger sequencing. The extracted nucleic acid is preferably sequenced using high-throughput next-generation sequencing (NGS). In principle, any paired-end sequencing method can be used to sequence the extracted DNA. In a preferred embodiment, direct target sequencing (DTS) is used, where sequences from the region of interest are enriched, if possible, using hybrid-capture probes or PCR primers designed such that the captured and sequenced fragments contain at least one sequence that distinguishes the target sequence from other captured sequences. In some embodiments, the paired-end reads obtained by DTS of one or many sites of interest include a first sequence read that includes a genome read and a second sequence read that includes a probe read associated with the site of interest in the genome of the subject. In some embodiments, the sequencing reads are 30-50 bp. In other embodiments, the sequencing reads are 100-200 bp in length. In a preferred embodiment, the sequence reads are about 40 bp. In some embodiments, DTS is used as described in U.S. Patent No. 9,092,401, which is incorporated by reference in its entirety.
[0046]
[0057] For example, hybrid-capture probes may be designed to anneal adjacent to a small number of bases that differ between different sites of interest ("diff bases"). If such discriminatory sequences are rare, multiple probes may be used to capture discriminable fragments, reducing the effect of bias inherent in the sequence of each particular probe.
[0047]
[0058] Nucleic acid sequences may be aligned to a reference genome to detect genetic variations. In a preferred embodiment, the subject is a human and the sequence reads are aligned to a human reference genome. For example, a sequence manipulation and alignment procedure ("pipeline") may begin with raw data from a genome analyzer, e.g., a Genome Analyzer IIx (GAIIx) or HiSeq sequencer (Illumina; San Diego, Calif.), to infer genotypes and calculate metrics from patient samples. Sequencing data from regions of interest may be sequenced according to the methods of the present invention using barcoded samples in a multiplexed (e.g., 12x) configuration per flow cell lane. The raw data from the sequencer may include base calls (BCL files) and various quality control and calibration metrics. The raw base calls and metrics may first be compiled into a QSEQ file, then filtered, merged, and demultiplexed (based on barcode sequences) into sample-specific FASTQ files. The FASTQ reads may be aligned to a reference genome, e.g., the HG19 genome, to create an initial BAM file. In some cases, each paired-end FASTQ file may be aligned to a reference genome. In other cases, each single-end FASTQ file may be aligned separately to a genome, allowing for "ambiguous alignments" and reporting of several top alignments for each read. In still other embodiments, the overall alignment process may include a single alignment of the forward and reverse paired-end NGS reads, and / or separate alignment or realignment of the forward and reverse single-end NGS reads (e.g., "ambiguous alignments"). The resulting BAM files can undergo several transformations to filter, clip, and refine alignments, and recalibrate quality metrics. The final BAM files can be used to estimate genotypes for known variants and discover new variants, resulting in a call set. The call set (VCF file) can then be filtered using various call metrics to produce a final set of reliable (e.g., about 80%, 85%, 90%, 95%, 99% or higher percentage confidence or greater than about 80%, 85%, 90%, 95%, 99% or higher percentage confidence) variant calls for each sample. Finally, various metrics can be calculated for each sample, lane, and batch, and the calls and metrics are loaded into a laboratory information management system (HMS) for visualization, review, and final report generation.The pipeline may be executed (in whole or in part) locally and / or using cloud computing such as in the Amazon cloud. A user may interact with the pipeline using any suitable communication mechanism. For example, interaction may be via Django administrative commands (Django Software Foundation, Lawrence, Kans.), shell scripts for executing each step of the pipeline, or an application programming interface written in a suitable programming language (e.g., an interface such as PHP, Ruby on Rails, Django, or Amazon EC2). An overview of the operation of this example pipeline is shown in Figures 10 and 11 of U.S. Pat. No. 9,092,401, which is incorporated herein by reference in its entirety.
[0048]
[0059] In some embodiments, alignment according to the present invention is performed using a computer program. One exemplary alignment program implementing the BWT approach is Burrows-Wheeler Aligner (BWA), available from the SourceForge website maintained by Geeknet (Fairfax, Va.). The quality of the alignments can be evaluated and / or compared by calculating an alignment score. For example, the quality of the alignments can be evaluated and / or compared by calculating an alignment score, as described by Heng Li (2013) "Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM" (arXiv:1303.3997v2[q-bio.GN]). The alignment score for each read or pair of reads can be used to identify a single top alignment or multiple top alignments for a collection of single-end or paired-end reads. In some cases, the aligner will only emit alignments that meet a minimum alignment score for a region of interest, e.g., the first, second, or subsequent regions of interest.
[0049]
[0060] Provided herein is a method for detecting genetic variation in a genome of a subject, wherein the genome comprises highly homologous regions of interest, and the detected genetic variation is within one or more of the highly homologous regions of interest.In some embodiments, the highly homologous regions have sequence identity of more than 70%, 71%, 72%, 73%, 74%, 75%, 76%, 77%, 78%, 79%, 80%, 81%, 82%, 83%, 84%, 85%, 86%, 87%, 88%, 89%, 90%, 91%, 92%, 93%, 94%, 95%, 96%, 97%, 98%, 99%, or 99%.In some cases, the method is effective for detecting genetic variation between two or more highly homologous regions in a genome.A highly homologous region may comprise any two or more regions that are highly similar. A homologous region may contain two or more genes that are highly similar. In some cases, a homologous region may contain one or more genes and one or more homologs of that gene. For example, a homolog may contain one or more pseudogenes. Genotyping such highly homologous regions using standard targeted NGS strategies, which use hybridization to capture and sequence short DNA fragments within each highly homologous region, is complicated by the fact that sequence reads cannot be unambiguously aligned to specific regions due to the relatively short read length and high homology between the regions. For example, PMS2 is usually included in HCS panels due to its association with Lynch syndrome [11-15]. Its nearby pseudogene, PMS2CL, complicates accurate NGS read alignment and identification of variants in exons 11 to 15 at the 3' end of PMS2 (Figure 1A): the coding sequence was previously reported to share 98% sequence identity with PMS2CL
[16] . Furthermore, sequence exchanges and gene conversions between the two regions are frequent enough that even a few non-identical bases in the reference genome (hg19) cannot be assigned with confidence to genes or pseudogenes [17, 18].Long-range PCR (LR-PCR) using gene-specific primers in exon 10 specifically amplifies PMS2 (Figure 1B), and variants in the terminal five exons of PMS2 can then be identified by Sanger sequencing [19-21] or NGS
[22] (Figure 1C). Although identification of copy number variants (CNVs) in PMS2 is possible from LR-PCR and Sanger sequencing, it is not straightforward, which motivated the parallel use of multiplex ligation-dependent probe amplification (MLPA) to detect large deletions and duplications [19-24].
[0050]
[0061] For highly homologous regions of the genome, e.g., PMS2, numerous testing strategies exist that can achieve high sensitivity and specificity ([18–20, 22, 25, 26]), but each requires quality control monitoring. For example, for the five final exons of PMS2, LR-PCR, MLPA, and hybrid-capture NGS on each screened sample have been previously published for small cohorts
[22] , but applying this combination to larger patient populations would be resource intensive and involve complex workflow logistics. Herman et al.
[26] recently presented a method for identifying CNVs (but not SNVs or indels) in the terminal exons of PMS2 or PMS2CL
[26] . This method identified samples for follow-up LR-PCR testing and, ultimately, localized CNVs to genes or pseudogenes. The authors noted a CNV false-positive rate of 6.8%, meaning that a significant portion of CNV-negative samples will undergo follow-up testing unnecessarily.
[0051]
[0062] High reflex rates after short-read NGS testing (e.g., >10%), although acceptable for patient reporting accuracy, can result in unmanageable logistical overhead at the testing laboratory. The reflex rate has two components, one biological component and one technical component, each with different sources and constraints. The biological component is used as the floor of the reflex rate. If the assay played a role as a marker for CNVs and had sufficient analytical specificity (i.e., zero false positives) and clinical accuracy (i.e., correct classification without VUS), then the reflex rate would nevertheless not be zero due to the presence of pathogenic variants in exons 12-15 of PMS2 and the corresponding PMS2CL region that requires disambiguation. This biological component therefore primarily reflects the cumulative population frequency of pathogenic variants across the ambiguity region. The technical component of the reflex rate, in contrast, results from insufficient analytical specificity and incomplete knowledge of variant pathogenicity. Although higher (99.7%) in Example 1, the analytical specificity for CNVs was 93.7% in Herman et al.
[26] , which meant that the technical component of the reflex rate in this study was at least 6.3% (highlighting the variable nature of the technical component). Furthermore, technical reflexes due to VUS in the workflow described herein were required in 4% of samples, an occupancy that is expected to decline with further screening of PMS2 and the resulting ability to reclassify VUS.
[0052]
[0063] Thus, a reflex method for detecting alterations between homologous regions in a genome is disclosed herein. The objective of the method is to have an initial testing phase of the workflow (i.e., upstream of the reflex) that is sensitive enough to maximize detection of PMS2 variants and specific enough to minimize reflex load. In one embodiment, the method applies hybrid-capture NGS to all samples and only LR-PCR / MLPA as the reflex assay. In some embodiments, the workflow described herein has high analytical accuracy (i.e., it is possible to detect sequence variants in specific regions), but requires reflex testing only for 10%, 9%, 8%, 7%, 6%, 5%, 4%, 3%, 2%, 1%, or less than 1% of the samples. In one embodiment, the workflow described herein has high analytical accuracy, but requires reflex testing only for about 8% of the samples. An exemplary embodiment of a method for detection of SNVs, indels, and CNVs in the five final exons of PMS2 is described in Example 1.
[0053]
[0064] In one embodiment, a method for detecting genetic variations in a genome of a subject, the genome comprising a first highly homologous region of interest and a second highly homologous region of interest, the method comprising: (a) obtaining sequence reads by paired-end sequencing from multiple sites of interest in the first and second regions of interest, the sequence reads comprising a first read and a second read obtained at each site of interest; (b) aligning the sequence reads to a reference genome, the first read and the second read being aligned separately to the reference genome, and an aligner generating multiple possible alignments for each of the first read and the second read; (c) identifying the first read and the second read that align to the first region of interest; (d) pairing the first read and the second read from the reads identified in step (c), thereby generating a top pair alignment; and (e) detecting genetic variations in the top pair alignment generated in step (d). In a preferred embodiment, the reads are aligned to a reference genome, the reference genome does not include masked or modified portions of the first homologous region or the second homologous region of interest, and the first homologous region and / or the second homologous region of interest are analyzed to detect genetic variations as described herein. The alignment in step (b) is referred to as "ambiguous alignment" because each single-end sequence read is aligned separately to the reference genome, and the alignment of multiple reads is identified in (c). An example of the implementation of the method according to the "ambiguous alignment" process is shown in Figure 9.
[0054]
[0065] In another embodiment, a method for detecting genetic variations in a genome of a subject, the genome comprising a first highly homologous region of interest and a second highly homologous region of interest, the method comprising the steps of: (a) obtaining sequence reads by paired-end sequencing from multiple sites of interest in the first and second regions of interest, the sequence reads comprising a first read and a second read obtained at each site of interest; and (b) aligning the first and second reads to a reference genome, the aligner generating a best possible paired-end alignment for each pair of the first and second reads to the first or second region of interest, and selecting a top aligner for the first or second region of interest. The method includes the steps of: (c) aligning the sequence reads to a reference genome, wherein the first and second reads are aligned to the reference genome separately, and the aligner generates multiple possible alignments for each of the first and second reads; (d) identifying the first and second reads that align to the first region of the object; (e) pairing the first and second reads from the reads identified in step (d), thereby generating a top pair alignment; and (f) detecting genetic variations in the top pair alignment generated in step (e). In a preferred embodiment, the reads are aligned to a reference genome, and the reference genome does not include a masked or modified portion of the first or second homologous region of the object, and the first and / or second homologous region of the object is analyzed to detect genetic variations as described herein. Thus, in some embodiments, a standard paired-end alignment is first performed to select reads that align to the region of interest, and typically only the paired-end reads with the top alignment scores are selected.The selected paired-end reads can then be partitioned and aligned separately to the reference genome to identify a number of top single-end alignments for each read (e.g., "ambiguous alignments").
[0055]
[0066] For each read, the multiple top single-end alignments generated by the aligner are individually paired to generate a top-paired alignment. For example, the top paired-end reads are partitioned into BAM files using, for example, samtool
[28] , and the BAM files are then parsed into a BAM file using, for example, Picard (Broadcast Each read is converted into two unaligned FASTQ files (each number of parsed read pairs into one of two files) using the Chromatography Institute, and each single-end FASTQ file is realigned separately against the reference genome, allowing for "ambiguous alignments," and the reporting of several top alignments for each read. Such top alignments can be used in the pairing step to identify the top pair alignment.
[0056]
[0067] The single-end reads selected from the "ambiguous alignment" may be used to generate a top paired-end alignment through a selection step. The single-end alignment may be used to generate a top paired-end alignment when: 1) both single-end reads have the same read name, 2) both single-end reads are mapped to a region spanning the region of interest used to identify single-end reads by "ambiguous alignment" as described above, and / or 3) both single-end reads align within a certain number of bases of each other. In a preferred embodiment, only reads that meet all of the pairing criteria (1)-(3) are paired. In some embodiments, the alignment of the first and second reads in the region of interest used to identify single-end reads by "ambiguous alignment" as described above is within about 100 bp, about 200 bp, about 300 bp, about 400 bp, about 500 bp, about 600 bp, about 700 bp, about 800 bp, about 900 bp, about 1000 bp, about 1200 bp, about 1400 bp, about 2000 bp, about 3000 bp, about 4000 bp, about 5000 bp, about 6000 bp, about 7000 bp, about 8000 bp, about 10000 bp, about 12000 bp, about 14000 bp, about 16000 bp, about 18000 bp, about 20000 bp, about 20000 bp, about 30000 bp, about 40000 bp, about 50000 bp, about 10000 bp, about 20000 bp, about 30000 bp, about 40000 bp, about 50000 bp, about 500 Reads are paired only if they are within about 100 bp, about 400 bp, about 500 bp, about 600 bp, about 700 bp, about 800 bp, about 900 bp, about 1000 bp, about 1100 bp, about 1200 bp, about 1300 bp, about 1400 bp, about 1500 bp, or more than 1500 bp. In some cases, if a large number of putative pairs meet the above conditions for a given read name, the pair with the highest alignment score is selected. In some cases, the top paired end alignment is selected as the one with the smallest template length. Reads that cannot form a suitable pair as described above are discarded. The resulting paired-end BAM file contains reads originating from both homologous regions of the target that are mapped to the region of the target used to identify the single-end reads by "ambiguous alignment". The top paired end alignment can be analyzed to identify or call variants in one or more homologous regions of the target.
[0057]
[0068] For example, for PMS2, the resulting single-end alignments can be used to generate paired-end alignments if the following criteria are met: 1) both single-end reads have the same read name, 2) both single-end reads map to a region spanning exons 12-15 of PMS2, 3) both single-end reads align within 1000bp of each other, 4) if multiple putative pairs meet the above conditions for a given read name, the pair with the highest alignment score is selected, and 5) reads that cannot form suitable pairs as above are discarded. The resulting paired-end BAM file contains reads originating from both PMS2 and PMS2CL reads mapped to the PMS2 sequence.
[0058]
[0069] In one embodiment, the genetic variation detected in the homologous sequence comprises one or more SNPs. In another embodiment, the genetic variation detected in the homologous sequence comprises one or more CNVs. In another embodiment, the genetic variation detected in the homologous sequence comprises one or more indels. In another embodiment, the genetic variation detected in the homologous sequence comprises one or more inversions. In another embodiment, the genetic variation detected in the homologous sequence comprises a combination of SNPs, indels, inversions, and / or CNVs.
[0059]
[0070] In one embodiment, for detecting genetic variations in the genome of a subject described herein, the genome comprises highly homologous regions including a first region and a second region of interest, and sequence reads are obtained from one or more exons in the first region and / or the second region of interest. The sequence reads can be obtained from one or more introns in the first region and / or the second region of interest. The sequence reads can be obtained from one or more exons and introns in the first region and / or the second region of interest. The sequence reads can be obtained from one or more exons and introns in the first region and / or the second region of interest, the introns being located near the exons. The introns located near the exons can be located within + / - 1 to 100 nt, e.g., + / - 20 nt, of the exons. The sequence reads can be obtained from one or more clinically addressable regions associated with the first region and / or the second region of interest. Such regions associated with the first and / or second regions of interest may include any region of the genome. For example, clinically addressable regions may include promoters, enhancers, and / or non-translated regions. In some cases, the first region of interest includes a gene and the second region of interest includes a pseudogene. In other cases, the first region of interest may include a pseudogene and the second region of interest includes a gene. The first region of interest may include two alleles. The second region of interest may include two alleles.
[0060]
[0071] In one embodiment, when a genetic variation is detected in a highly homologous region of interest in a genome of a subject by the methods described herein, a portion of the genome of the subject is amplified by long-range PCR and assayed by multiplex ligation-dependent probe amplification (MLPA). In another embodiment, when a genetic variation is detected in a highly homologous region of interest in a genome of a subject by the methods described herein, a portion of the first region of interest is amplified by long-range PCR and the product or a portion thereof is sequenced by Sanger sequencing. In another embodiment, when a genetic variation is detected in a highly homologous region of interest in a genome of a subject by the methods described herein, a portion of the first region of interest is amplified by long-range PCR and the product or a portion thereof is sequenced by NGS. In another embodiment, when a genetic variation is detected in a highly homologous region of interest in a genome of a subject by the methods described herein, the genomic DNA of the subject is assayed by multiplex ligation-dependent probe amplification (MLPA).
[0061]
[0072] In one embodiment, the gene is PMS2 and the pseudogene is PMS2CL or one of several other pseudogenes for PMS2. The pseudogenes for exons 9 and 11-15 of PMS2 may be selected from, but are not limited to, PMS2CL. The pseudogenes for all of PMS2, particularly exons 1-5 of PMS2, may be selected from, but are not limited to, 15 or more / fewer pseudogenes. In an embodiment, the presence of altered copy numbers and / or inversions that change the orientation of the gene and pseudogene (e.g., fusing a portion of the pseudogene with the gene, thus compromising the function of the gene) may indicate that the subject has an increased risk for the disease Lynch syndrome.
[0062]
[0073] In one embodiment, the multiple sites of interest in the highly homologous region from which paired-end reads are obtained are within an exon of PMS2 and an exon of another part of the genome of the subject. In another embodiment, the multiple sites of interest are within an exon of PMS2 and an exon of PMS2CL. In another embodiment, the multiple sites of interest are within exons 11, 12, 13, 14, and / or 15 of PMS2 and exons 2, 3, 4, 5, and / or 6 of PMS2CL.
[0063]
[0074] In one embodiment, the gene is SMN1 and the pseudogene is SMN2. In an embodiment, the presence of an altered copy number of SMN1 indicates that the subject may be a carrier for the disease spinal muscular atrophy (SMA).
[0064]
[0075] In another embodiment, the gene is CYP21A2 and the pseudogene is CYP21A1P. In an embodiment, the presence of an altered copy number of CYP21A2 indicates that the subject may be a carrier for the disease congenital adrenal hyperplasia (CAH).
[0065]
[0076] In an embodiment, the gene is HBA1 and the homolog is HBA2 (or vice versa). In an embodiment, the presence of an altered copy number of either HBA1 or HBA2 indicates that the subject may be a carrier for the disease alpha thalassemia.
[0066]
[0077] In a further embodiment, the gene is GBA and the pseudogene is GBAP. In an embodiment, the presence of an altered copy number of GBA indicates that the subject may be a carrier for the disease Gaucher's disease.
[0067]
[0078] In an embodiment, the gene is CHEK2, which has several pseudogenes. As of December '14, there were seven pseudogenes. The pseudogene may be selected from, but is not limited to, the CHEK2 pseudogenes listed in curated databases. In embodiments, the presence of a mutation resulting from recombination with the pseudogene, e.g., a frameshift mutation originating from the pseudogene, may indicate that the subject has an increased risk for the disease breast cancer, among other diseases. It is well known in the art that only one of the seven pseudogenes has been named, and that the risk is primarily associated with one mutation, namely 1100delC. However, other mutations also contribute to the risk of the disease. The patient is at risk for Li-Fraumeni syndrome and other hereditary cancers.
[0068]
[0079] In embodiments, the gene is SDHA and the pseudogene is any one of the pseudogenes thereof, for example, SDHAP1, SDHAP2, SDHAP3.
[0069] III. Variant calling
[0080] In some embodiments, variants are detected with a computer-implemented caller algorithm. In principle, any variant caller can be utilized, for example, to detect SNPs, indels, inversions, and CNVs. In some cases, a caller is used that is capable of detecting / resolving breakpoints when genetic variations, such as deletions, are detected. For example, the caller can be selected from the callers described in Tattini, L. et al., Front Bioeng Biotechnol. 2015; 3: 92. In some cases, variants are identified based on an expected ploidy of 0-7, or 0-8. In some cases, variants are identified based on an expected ploidy of 2. In other cases, variants are identified based on an expected ploidy of 6. In other cases, variants are identified based on an expected ploidy of 4. For example, SNVs and indels can be identified using GATK 4.0 HaplotypeCaller
[29] with the sample-ploidy option set to 4 (e.g., for the exon 12-15 region of tetraploid PMS2). In other cases, SNVs and short indels can be identified using GATK 1.6
[30] and FreeBayes
[31] with the sample-ploidy option set to 2 (e.g., for the exon 11 region of diploid PMS2). For diploid SNV calling in LR-PCR data, GATK 1.6 can be used similarly.
[0070]
[0081] In a preferred embodiment, a Hidden Markov Model (HMM) caller is used to determine copy number. The preferred caller used to determine copy number is the HMM caller described in U.S. Provisional Patent Application No. 62 / 681,517, which is incorporated herein by reference in its entirety. In some embodiments, the preferred HMM caller is set to an expected ploidy of 2. In other embodiments, the preferred HMM caller is set to an expected ploidy of 4. In other embodiments, the preferred HMM caller is set to an expected ploidy of 6.
[0071]
[0082] In some embodiments, methods for evaluating sample-specific performance of copy number variant models, methods for determining copy numbers of interrogated segments within a region of interest, and methods for determining copy number variant aberrations within a region of interest are utilized as described in U.S. Provisional Patent Application No. 62 / 681,517, which is incorporated by reference in its entirety herein.
[0072]
[0083] In some embodiments, a method for evaluating sample-specific performance of a copy number variant caller comprising a copy number variant model includes the steps of: parameterizing the copy number variant model based on an actual number of sequencing reads mapped to a segment within a region of interest from a test sample, determining one or more copy number variant model parameters; and generating a plurality of synthetic copy number variants. wherein each synthetic copy number variant comprises one or more synthetic copy numbers of a segment, each synthetic copy number represented by a synthetic sequencing read number based on an actual number of sequencing reads for a corresponding segment from the test sample; calling copy numbers of the one or more segments for the synthetic copy number variant using a copy number variant model and one or more determined copy number variant model parameters; determining a sample-specific performance statistic for the copy number variant caller based on a difference between the called copy numbers and the synthetic copy numbers in the synthetic copy number variants; and evaluating the sample-specific performance of the copy number variant caller based on the sample-specific performance statistic.
[0073]
[0084] In some embodiments of the method for evaluating the sample-specific performance of copy number variant caller, the number of synthetic sequencing reads for one or more segments is obtained by increasing, decreasing or maintaining the number of actual sequencing reads for the corresponding segments from test samples in proportion to a predetermined number of copies of one or more segments.In some embodiments, the predetermined number of copies is an integer number of copies.In some embodiments, the predetermined number of copies is a non-integer number of copies.
[0074]
[0085] In some embodiments of the method of evaluating the sample-specific performance of a copy number variant caller, the synthetic sequencing read number is obtained by sampling a binomial distribution with probability of success equal to m / x and a number of tests equal to the actual sequencing read number for the corresponding segment from the test sample, where m is the synthetic copy number of the segment in the synthetic copy number variant and x is the hypothetical copy number of the corresponding segment from the test sample.
[0075]
[0086] In some embodiments of the method for evaluating the sample-specific performance of a copy number variant caller, the synthetic sequencing read number is obtained by sampling the number of sequencing reads as a negative binomial distribution with a probability of success equal to m / x and a number of successes equal to the actual number of sequencing reads in the corresponding segment from the test sample, where m is the synthetic copy number of the segment in the synthetic copy number variant and x is the hypothetical copy number of the corresponding segment from the test sample, and adding the sampled sequencing read number to the actual sequencing read number for the corresponding segment from the test sample. In some embodiments, the synthetic sequencing read number is obtained by sampling the number of sequencing reads as an expectation of the negative binomial distribution.
[0076]
[0087] In some embodiments of the method for evaluating sample-specific performance of a copy number variant caller, the copy number variant model is a hidden Markov model. In some embodiments, the hidden Markov model includes: (i) one or more hidden states including copy numbers corresponding to the interrogated segment or multiple subsegments within the interrogated segment, (ii) an observed state including the number of actual or synthetic sequencing reads for the interrogated segment, and (iii) a copy number likelihood model based on the predicted number of actual or synthetic sequencing reads for the interrogated segment. In some embodiments, the method includes determining the copy number likelihood model. In some embodiments, parameterizing the hidden Markov model includes adjusting the copy number likelihood model to fit the actual number of sequencing reads mapped to the interrogated segment from the test sample. In some embodiments, the copy number likelihood model includes a distribution for two or more copy number states. In some embodiments, the copy number likelihood model includes a negative binomial distribution, where the negative binomial distribution is not a Poisson distribution. In some embodiments, the predicted number of actual or synthetic sequencing reads is based on the segmental states corresponding to the interrogated segment across multiple samples. Based on the average number of sequencing reads mapped in the segment corresponding to the interrogated segment across the multiple samples, and the average number of sequencing reads mapped across the segments in the test sample, the average number of sequencing reads mapped in the segment corresponding to the interrogated segment across the multiple samples or the average number of sequencing reads mapped across the multiple segments in the test sample is a normalized average. In some embodiments, the copy number likelihood model is adjusted to take into account the presence of GC content bias. In some embodiments, the hidden Markov model includes a transition probability of the copy number of the interrogated segment for a given copy number of a spatially adjacent segment. In some embodiments, the hidden Markov model includes a plurality of transition probabilities of the copy number of the subsegment in the interrogated segment for a given copy number of a spatially adjacent subsegment. In some embodiments, the transition probability takes into account the average length of the copy number variant. In some embodiments, the transition probability takes into account the prior probability of the copy number variant in the interrogated segment or in the spatially adjacent segments. In some embodiments, the average length of the copy number variant or the probability of the copy number variant in the interrogated segment is determined based on observations in a human population.
[0077]
[0088] In some embodiments of the method for evaluating the sample-specific performance of a copy number variant caller, parameterizing the copy number variant model comprises taking into account one or more false capture probes. In some embodiments, taking into account one or more false capture probes comprises weighting one or more observation states in a plurality of observation states including a false capture probe indicator. In some embodiments, the false capture probe indicator is determined using a Bernoulli process. In some embodiments, taking into account one or more of the false capture probes comprises using expectation maximization. In some embodiments, if a capture probe is determined to be false, the sequencing reads from that capture probe are ignored in the copy number variant model.
[0078]
[0089] In some embodiments of the method of evaluating the sample-specific performance of a copy number variant caller, parameterizing the copy number variant model includes taking into account noise in the number of mapped sequencing reads.
[0079]
[0090] In some embodiments of the methods of assessing sample-specific performance of a copy number variant caller, the copy number variant model is parameterized using the analytical gradient of a first derivative and a Hessian of a second derivative of one or more copy number variant model parameters.
[0080]
[0091] In some embodiments of the methods of evaluating the sample-specific performance of a copy number variant caller, the copy number variant model is parameterized by solving a trust region Newton conjugate gradient algorithm.
[0081]
[0092] In some embodiments of the methods of evaluating the sample-specific performance of a copy number variant caller, the copy number variant model is iteratively parameterized using expectation maximization.
[0082]
[0093] In some embodiments of a method for evaluating sample-specific performance of a copy number variant caller, the method includes mapping actual sequencing reads from a test sample to segments within a region of interest and determining the actual number of sequencing reads mapped to the segments.
[0083]
[0094] Some embodiments of methods for assessing sample-specific performance of a copy number variant caller In, a test sample is enriched using one or more direct targeted sequencing capture probes.
[0084]
[0095] In some embodiments of the method of evaluating the sample-specific performance of a copy number variant caller, the method includes calling the copy number of one or more segments for a test sample.
[0085]
[0096] In some embodiments of the methods of evaluating the sample-specific performance of a copy number variant caller, the segments comprise spatially adjacent segments.
[0086]
[0097] In some embodiments of the methods of assessing the sample-specific performance of a copy number variant caller, the sample-specific performance statistic is limit of detection, sensitivity, specificity, accuracy, recall, precision, positive predictive value, or negative predictive value.
[0087]
[0098] In some embodiments of the methods of assessing the sample-specific performance of a copy number variant caller, the sample-specific performance statistic is sensitivity or precision.
[0088]
[0099] In some embodiments of the method of evaluating the sample-specific performance of a copy number variant caller, the method includes rejecting the test sample if the sample-specific performance of the copy number variant model is below a desired performance threshold.
[0089]
[0100] A method for determining the copy number of an interrogated segment within a region of interest, comprising: (a) mapping a plurality of sequencing reads generated from a test sequencing library to the interrogated segment, the test sequencing library being enriched using one or more direct targeted sequencing capture probes; (b) determining a number of sequencing reads mapped to the interrogated segment; (c) determining a copy number likelihood model based on a predicted number of sequencing reads mapped to the interrogated segment; and (d) generating one or more hidden state likelihood models comprising (i) copy numbers corresponding to the interrogated segment or a plurality of subsegments within the interrogated segment. Also described herein is a method comprising: (i) constructing a hidden Markov model including an observation state, the observation state including the number of sequencing reads mapped to the interrogated segment, and (ii) a copy number likelihood model; (e) parameterizing the hidden Markov model by adjusting the copy number likelihood model to fit a predetermined number of sequencing reads mapped to the interrogated segment, the hidden Markov model being parameterized using an analytical gradient of the first derivative and a Hessian matrix of a second derivative of one or more parameters in the copy number likelihood model; and (f) determining the most likely copy number of the interrogated segment based on the parameterized hidden Markov model.
[0090]
[0101] A method for determining the copy number of an interrogated segment within a region of interest, comprising: (a) mapping a plurality of sequencing reads generated from a test sequencing library to a plurality of spatially adjacent segments, the plurality of spatially adjacent segments including the interrogated segment, the test sequencing library being enriched using a plurality of spatially adjacent direct target sequencing capture probes; (b) determining a number of sequencing reads mapped to each spatially adjacent segment; (c) determining a copy number likelihood model for each spatially adjacent segment based on a predicted number of sequencing reads mapped in the spatially adjacent segment; and (d) (i) determining a copy number likelihood model for each of the spatially adjacent segments or within each of the spatially adjacent segments. Further described herein is a method comprising: constructing a hidden Markov model including (i) a plurality of hidden states including copy numbers for each of a plurality of subsegments to be investigated; (ii) a plurality of observed states including a number of sequencing reads mapped to each spatially proximate segment; and (iii) a copy number likelihood model for each spatially proximate segment; (e) parameterizing the hidden Markov model including adjusting each copy number likelihood model to fit a predetermined number of sequencing reads mapped to each spatially proximate segment, wherein the hidden Markov model is parameterized using an analytical gradient of a first derivative and a Hessian matrix of a second derivative of one or more parameters in the copy number likelihood model; and (f) determining a most likely copy number for the interrogated segment based on the parameterized hidden Markov model.
[0091]
[0102] A method for determining copy number variant aberrations in a region of interest, comprising: (a) mapping a plurality of sequencing reads generated from a test sequencing library to an interrogated segment in the region of interest, the test sequencing library being enriched using one or more direct targeted sequencing capture probes; (b) determining a number of sequencing reads mapped to the interrogated segment; (c) determining a copy number likelihood model based on an expected number of sequencing reads mapped to the interrogated segment; and (d) determining a copy number likelihood model based on an expected number of sequencing reads mapped to the interrogated segment. Also described herein are methods that include constructing a hidden Markov model including an observation state including a number of sequencing reads mapped to the interrogated segment, and (iii) a copy number likelihood model; (e) parameterizing the hidden Markov model by adjusting the copy number likelihood model to fit a predetermined number of sequencing reads mapped to the interrogated segment, wherein the hidden Markov model is parameterized using an analytical gradient of a first derivative and a Hessian matrix of a second derivative of one or more parameters in the copy number likelihood model; (f) determining a most likely copy number for the interrogated segment based on the parameterized hidden Markov model; and (g) determining a copy number variant aberration based on the most likely copy number for the interrogated segment.
[0092]
[0103] A method for determining copy number variant abnormalities in a region of interest, comprising: (a) mapping a plurality of sequencing reads generated from a test sequencing library to a plurality of spatially adjacent segments, the plurality of spatially adjacent segments including an interrogated segment, the test sequencing library being enriched using a plurality of spatially adjacent direct target sequencing capture probes; (b) determining a number of sequencing reads mapped to each spatially adjacent segment; (c) determining a copy number likelihood model for each spatially adjacent segment based on a predicted number of sequencing reads mapped in the spatially adjacent segment; and (d) determining a copy number likelihood model for each spatially adjacent segment based on a predicted number of sequencing reads mapped in the spatially adjacent segment. constructing a hidden Markov model including (i) a plurality of hidden states including copy numbers for each of a plurality of subsegments within each of the segments or each of the spatially adjacent segments, (ii) a plurality of observed states including a number of sequencing reads mapped to each of the spatially adjacent segments, and (iii) a copy number likelihood model for each of the spatially adjacent segments; and (e) parameterizing the hidden Markov model including adjusting each copy number likelihood model to fit a predetermined number of sequencing reads mapped to each of the spatially adjacent segments, wherein the hidden Markov model is parameterized using an analytical gradient of a first derivative and a Hessian of a second derivative of one or more parameters in the copy number likelihood model. Further described herein is a method comprising the steps of: (a) determining a most likely copy number of the interrogated segment based on a parameterized hidden Markov model; and (b) determining copy number variant abnormalities based on the most likely copy number of the interrogated segment.
[0093]
[0104] A method for determining the copy number of an interrogated segment within a region of interest, comprising: (a) mapping a plurality of sequencing reads generated from a test sequencing library to the interrogated segment, the test sequencing library being enriched using one or more capture probes; (b) determining a number of sequencing reads mapped to the interrogated segment; (c) determining a copy number likelihood model based on a predicted number of sequencing reads mapped to the interrogated segment; and (d) determining a copy number likelihood model based on a predicted number of sequencing reads mapped to the interrogated segment. Also described herein are methods comprising: (i) constructing a hidden Markov model comprising an observation state including a number of sequencing reads mapped to the interrogated segment, and (ii) a copy number likelihood model; (e) parameterizing the hidden Markov model by adjusting the copy number likelihood model to fit a predetermined number of sequencing reads mapped to the interrogated segment and taking into account one or more false capture probes, wherein the hidden Markov model is parameterized using an analytical gradient of a first derivative and a Hessian matrix of a second derivative of one or more parameters in the copy number likelihood model; and (f) determining the most likely copy number of the interrogated segment based on the parameterized hidden Markov model.
[0094]
[0105] A method for determining the copy number of an interrogated segment within a region of interest, comprising: (a) mapping a plurality of sequencing reads generated from a test sequencing library to a plurality of spatially adjacent segments, the plurality of spatially adjacent segments including the interrogated segment, the test sequencing library being enriched using a plurality of spatially adjacent direct target sequencing capture probes; (b) determining a number of sequencing reads mapped to each spatially adjacent segment; (c) determining a copy number likelihood model for each spatially adjacent segment based on a predicted number of sequencing reads mapped in the spatially adjacent segment; and (d) (i) determining a copy number likelihood model for each of the spatially adjacent segments or each of a plurality of subsegments within each of the spatially adjacent segments. Further described herein is a method comprising: constructing a hidden Markov model including (i) a plurality of hidden states including a number of sequences mapped to each spatially adjacent segment, (ii) a plurality of observed states including a number of sequencing reads mapped to each spatially adjacent segment, and (iii) a copy number likelihood model for each spatially adjacent segment; (e) parameterizing the hidden Markov model including adjusting each copy number likelihood model to fit a predetermined number of sequencing reads mapped to each spatially adjacent segment and taking into account one or more false capture probes, wherein the hidden Markov model is parameterized using an analytical gradient of a first derivative and a Hessian matrix of a second derivative of one or more parameters in the copy number likelihood model; and (f) determining the most likely copy number of the interrogated segment based on the parameterized hidden Markov model.
[0095]
[0106] In some embodiments of the above methods, one or more parameters of the copy number likelihood model are determined based on the variance (d i), the average number of mapped sequencing reads for a segment (μ i ), some mapping to segments in the test sequencing library Sequencing read variance (d j ), or the average number of mapped sequencing reads to a segment in the test sequencing library (μ j ).
[0096]
[0107] In some embodiments of the above method, the method further comprises a step of determining the most probable copy number of a section within the region of interest, the section comprising a plurality of spatially adjacent segments including the interrogated segment.
[0097]
[0108] In some embodiments of the above methods, the copy number likelihood model comprises distributions for two or more copy number states.
[0098]
[0109] In some embodiments of the above methods, the copy number likelihood model comprises a negative binomial distribution that is not a Poisson distribution.
[0099]
[0110] In some embodiments of the above method, the predicted number of sequencing reads is based on a normalized average of the average number of sequencing reads mapped in the corresponding segment across the multiple sequencing libraries and the average number of mapped sequencing reads across the multiple segments of the target in the test sequencing library.
[0100]
[0111] In some embodiments of the above methods, the copy number likelihood model is adjusted to take into account the presence of GC content bias. In some embodiments, the adjustment varies depending on the GC content of the capture probe corresponding to the interrogated segment or the GC content of the interrogated segment.
[0101]
[0112] In some embodiments of the above method, the hidden Markov model includes a transition probability of the copy number of the investigated segment for a given copy number of a spatially adjacent segment. In some embodiments, the transition probability takes into account the average length of the copy number variant. In some embodiments, the transition probability takes into account the previous probability of the copy number variant in the investigated segment or the spatially adjacent segment. In some embodiments, the average length of the copy number variant or the probability of the copy number variant in the investigated segment is determined based on observations in a human population.
[0102]
[0113] In some embodiments of the above method, the hidden Markov model includes a plurality of transition probabilities of copy numbers of subsegments in a plurality of subsegments in the interrogated segment for a given copy number of a spatially adjacent subsegment. In some embodiments, the transition probabilities take into account the average length of copy number variants. In some embodiments, the transition probabilities take into account the prior probability of copy number variants in the interrogated segment or spatially adjacent segments. In some embodiments, the average length of copy number variants or the probability of copy number variants in the interrogated segment are determined based on observations in a human population.
[0103]
[0114] In some embodiments of the above method, parameterizing the hidden Markov model includes taking into account one or more false capture probes. In some embodiments, taking into account one or more false capture probes includes weighting one or more observation states in a plurality of observation states that include a false capture probe indicator. In some embodiments, the false capture probe indicator is determined using a Bernoulli process. In some embodiments, taking into account one or more of the false capture probes includes using expectation maximization. In some embodiments, if a capture probe is determined to be false, the likelihood information from that capture probe is ignored in the copy number likelihood model.
[0104]
[0115] In some embodiments of the above method, parameterizing the hidden Markov model includes taking into account noise in the number of mapped sequencing reads.
[0105]
[0116] In some embodiments of the above method, taking into account the noise of the number of mapped sequencing reads comprises adjusting a copy number likelihood model. In some embodiments, adjusting the copy number likelihood model to take into account the noise comprises an expectation maximization step. In some embodiments, the expectation maximization step comprises weighting the level of noise of the number of mapped sequencing reads from the test sequencing library. In some embodiments, the most possible copy number of the investigated segment is not called if the noise of the number of mapped sequencing reads exceeds a predetermined threshold.
[0106]
[0117] In some embodiments of the above methods, sequencing reads from overlapping capture probes are merged.
[0107]
[0118] Some embodiments of the above method use a Viterbi algorithm, a quasi-Newton solver, or a Markov Chain Monte Carlo method to determine the most likely copy number of the interrogated segment.
[0108]
[0119] In some embodiments of the above methods, the method further comprises determining the reliability of the most possible copy number of the segment.
[0109]
[0120] In some embodiments of the above methods, one or more parameters of the copy number likelihood model are determined based on the variance (d i ), the average number of mapped sequencing reads for a segment (μ i ), the distribution of several mapped sequencing reads to segments in the test sequencing library (d j), or the average number of mapped sequencing reads to a segment in the test sequencing library (μ j ).
[0110]
[0121] In some embodiments of the above method, the analytical gradient of the first derivative and the analytical Hessian of the second derivative of the one or more parameters in the copy number likelihood model are solved using a trust-region Newton conjugate gradient algorithm.
[0111]
[0122] Also described herein is a computer system that includes a computer readable medium containing instructions for performing any one of the above methods.
[0112] IV. EXEMPLARY ARCHITECTURE AND PROCESSING ENVIRONMENT
[0123] In a preferred embodiment, some of the methods described herein are computer-implemented. Exemplary environments and systems in which certain aspects and examples of the systems and processes described herein may operate. As shown in FIG. 10, in some examples, the system may be implemented according to a client-server model. The system may include a client-side portion that runs on a user device 102 and a server-side portion that runs on a server system 110. The user device 102 may include any electronic device, such as a desktop computer, a laptop computer, a tablet computer, a PDA, a mobile phone (e.g., a smart phone), etc.
[0113]
[0124] The user device 102 may be connected to a network, including the Internet, an intranet, or any other wired or wireless public or private network. The client-side portion of the exemplary system on the user device 102 may provide client-side functionality, e.g., user-facing input and output processing and communication with the server system 110. The server system 110 may provide server-side functionality for any number of clients residing on each user device 102. Additionally, the server system 110 may include one or more caller servers 114, which may include a client-facing I / O interface 122, one or more processing modules 118, a data and model storage 120, and an I / O interface 116 to external services. The client-facing I / O interface 122 may facilitate client-facing input and output processing for the caller server 114. The one or more processing modules 118 may include various problem and candidate scoring models as described herein. In some examples, the caller server 114 may communicate with external services 124, e.g., text databases, subscription services, government records services, etc., through the network 108 for task completion or information acquisition. An I / O interface 116 for external services can facilitate such communication.
[0114]
[0125] The server system 110 can be implemented on one or more standalone data processing devices or distributed computer networks. In some examples, the server system 110 can use various virtual devices and / or services of a third party service provider (e.g., a third party cloud service provider) to provide the underlying computational and / or infrastructure resources of the server system 110.
[0115]
[0126] Although the functionality of the caller server 114 is illustrated in FIG. 10 as including both client-side and server-side portions, in some examples, certain functionality described herein (e.g., relating to user interface features and graphical elements) may be implemented as a standalone application installed on the user device. Furthermore, the division of functionality between the client and server portions of the system may vary in different examples. For example, in some examples, a client executing on the user device 102 may be a thin client that provides only user-facing input and output processing functions and delegates all other functionality of the system to a back-end server.
[0116]
[0127] It should be noted that the server system 110 and the client 102 may further include any one of various types of computing devices having, for example, a processing unit, memory (which may include logic or software for performing some or all of the functions described herein), and a communication interface, as well as other conventional computing components (e.g., input devices such as keyboard / touch screen, and output devices such as display). Furthermore, one or both of the server system 110 and the client 102 are generally programmed to include logic (e.g., http web server logic) or to format data accessed from local or remote databases or other data and content sources. To this end, the server system 110 may utilize various web data interface techniques, such as the Common Gateway Interface (CGI) protocol and associated applications (or "scripts"), Java "servlets", i.e., Java applications executing on the server system 110, to present information and receive input from the client 102. Although server system 110 is referred to in the singular herein, in reality it may comprise multiple devices (wired and / or wireless) for performing some or all of the functionality described herein. The server system 110 may include multiple computers, devices, databases, associated back-end devices, etc. that communicate and cooperate (over a line). The server system 110 may further include or communicate with an account server (e.g., an email server), a mobile server, a media server, etc.
[0117]
[0128] Additionally, while the exemplary methods and systems described herein describe the use of separate server and database systems to perform various functions, it should be noted that other embodiments may be implemented by storing software or programming that operates to cause the described functions on a single device or any combination of multiple devices as a matter of design choice, so long as the described functionality is performed. Similarly, the described database systems may be implemented as a single database, a distributed database, a collection of distributed databases, databases with redundant online or offline backups or other redundancies, and the like, and may include a distributed database or storage network and associated processing intelligence. Although not shown in the figures, the server system 110 (and other servers and services described herein) generally includes art-recognized components as typically found in a server system, including, but not limited to, a processor, RAM, ROM, clock, hardware drivers, associated storage, and the like (see, for example, FIG. 11 discussed below). Additionally, the described functions and logic may be included in software, hardware, firmware, or a combination thereof.
[0118]
[0129] FIG. 11 illustrates an exemplary computing system 1400 configured to perform any one of the above processes, including various calling and scoring models. In this context, the computing system 1400 may include, for example, a processor, memory, storage, and input / output devices (e.g., a monitor, keyboard, disk drives, an Internet connection, etc.). However, the computing system 1400 may include circuitry or other dedicated hardware for performing some or all aspects of the process. In some operating environments, the computing system 1400 may be configured as a system including one or more units, each configured to perform some aspects of the process, either in software, hardware, or some combination thereof.
[0119]
[0130] 11 illustrates a computing system 1400 with several components that may be used to implement the above process. The main system 1402 includes a motherboard 1404 having an input / output ("I / O") section 1406, one or more central processing units ("CPUs") 1408, and a memory section 1410 that may have associated therewith a flash memory card 1412. The I / O section 1406 is connected to a display 1424, a keyboard 1414, a disk storage unit 1416, and a media drive unit 1418. The media drive unit 1418 may read / write to a computer readable medium 1420 that may store programs 1422 and / or data.
[0120]
[0131] At least some values based on the results of the above processes can be stored for later use. Furthermore, a non-transitory computer-readable storage medium can be used to store (e.g., tangibly embody) one or more computer programs for performing any one of the above processes by a computer. The computer programs can be written, for example, in a general-purpose programming language (e.g., Pascal, C, C++, Python, Java) or some dedicated application-specific language.
[0121]
[0132] Various exemplary embodiments are described herein. Reference is made to these examples in a non-limiting sense. They are provided to illustrate the more widely applicable aspects of the disclosed technology. Various changes may be made and equivalents may be substituted without departing from the true spirit and scope of the various embodiments. Furthermore, many modifications may be made to adapt a particular situation, material, composition of matter, process, process act or step to the goal, spirit or scope of the various embodiments. Furthermore, as will be recognized by those skilled in the art, each of the individual variations described and illustrated herein has individual components and features that can be easily separated from or combined with the features of any of the other several embodiments without departing from the scope or spirit of the various embodiments. All such modifications are intended to be within the scope of the claims associated with this disclosure. EXAMPLES
[0122]
[0133] The present invention is further detailed in the following examples, which are not intended to limit the scope of the claimed invention in any way. The attached figures are meant to be considered as an integral part of the specification and description of the present invention. All references cited are specifically incorporated herein by reference for all that is stated therein. The following examples are provided to illustrate, but not to limit, the claimed invention.
[0123] Example 1 Detection of clinically actionable variants in the 3' exon of PMS2
[0134] This example shows a strategy for the detection of SNVs, indels, and CNVs in the 3' exon of PMS2. This study was reviewed and designated as exempt by the Western Institutional Review Board and complied with the Health Insurance Portability and Accountability Act (HIPAA).
[0124] material and method Research Samples:
[0135] Supplementary Table S1 indicates which sample sets were used for specific assays and analyses. Cell line DNA was purchased from Coriell Cell Repositories (Camden, NJ) (Supplementary Table S2). Patient sample DNA was extracted from de-identified blood or saliva samples. DNA samples with known positives were donations from Invitae Corporation.
[0125] LR-PCR:
[0136] DNA was extracted and further purified by incubation with 1x SPRI beads, followed by washing with 80% ethanol and elution in TE (10 mM Tris-HCl, 1 mM EDTA, pH 8.0). Approximately 300 ng of eluted DNA served as template in separate gene- and pseudogene-specific LR-PCR reactions with the following final concentrations: 1x LongAmp Taq Reaction Buffer (New England Biolabs, NEB), 0.3 mM dNTPs, 1 μM gene- or pseudogene-specific forward primer, 1 μM common reverse primer LRPCR_Unv_R (all primer sequences in the accompanying Table S3), 0.25% formamide, and 5 units of LongAmp Hot Start Taq Reactions containing the gene-specific forward primer PMS2_LRPCR_F yielded an approximately 17 kb amplicon spanning exons 11–15 of PMS2 (forward primer targeting exon 10), whereas use of the pseudogene-specific forward primer PMS2CL_F yielded an approximately 17 kb amplicon spanning exons 11–15 of PMS2 (forward primer targeting exon 10). Approximately 18 kb was amplified from the 5' end of pGFR2 (spanning the region upstream of PMS2CL from pGFR2.6). Thermal cycling included an initial denaturation at 94°C for 5 min, followed by 30 cycles of 94°C for 30 s and 65°C for 18.5 min. Final extension was at 65°C for 18.5 min, followed by a hold at 4°C. The quality of the LR-PCR amplicons was assessed using 0.5% agarose gel electrophoresis and quantified by a broad range Qubit assay kit (Thermo Fisher).
[0126]
[0137] Two different library prep strategies were used to prepare LR-PCR amplicons for NGS. First, for application to patient samples, LR-PCR amplicons were fragmented by adding 2 μL of NEBNext dsDNA Fragmentase and NEBNext dsDNA Fragmentase Reaction Buffer v2 (1× final, NEB) to the remaining LR-PCR reaction volume, then incubated at 37° C. for 25 min. The reaction was stopped by adding 100 mM EDTA and purified using 1.5× SPRI beads, followed by washing with 80% ethanol and elution in TE. The quality of fragmentation was assessed by Bioanalyzer (Agilent) using the High Sensitivity DNA kit. NGS library prep included end repair, A-tailing, and adapter ligation. Samples were PCR amplified for 8–10 cycles using the KAPA HiFi HotStart PCR Kit (Kapa Biosystems) containing barcoded primers with the following thermal cycling: initial denaturation at 95°C for 5 min, followed by cycles of 98°C for 20 s, 60°C for 30 s, and 72°C for 30 s. Final extension was at 72°C for 5 min, followed by a hold at 4°C. Library quality was assessed by Bioanalyzer using the High Sensitivity DNA kit, and concentration was measured by absorbance with a microplate reader (Tecan Infinite M200 PRO).
[0127]
[0138] The second approach to prepare LR-PCR amplicons for NGS was applied to samples of 155 cell lines and involved fragmenting and inserting adapters into LR-PCR amplicons by tagmentation. Two double-stranded adapters were generated by annealing single-stranded oligonucleotides: one double-stranded adapter had Unv_Tn5_oligo (all primer sequences in Table S3) annealed to Oligo A; the other double-stranded adapter had Unv_Tn5_oligo annealed to Oligo B. Two separate annealing mixes contained 25 μM of each oligonucleotide in double stranded and 1× annealing buffer (10 mM Tris-HCl, 50 mM NaCl, 1 mM EDTA, pH 8.0). The reactions were denatured at 95°C for 2 min and incubated at 80°C for 60 min, decreasing the temperature by 1 degree per minute until reaching 20°C, then kept at 4°C. Adapters were loaded into the Tn5 enzyme during a 30 min incubation at 37 °C using 0.15 units of Robust Tn5 Transposase (kit from Creative Biogene), 1.25 μM of each adapter, and 1 × TPS buffer. LR-PCR amplicons were subjected to tagmentation with the Tn5 adapter constructs. Using 0.5 μL of loaded Tn5 and 1–2 ng of DNA from each LR-PCR reaction, the tagmentation reaction was allowed to occur for 10 min at 56 °C in 1 × LM buffer. After incubation, SDS (0.02% final) was added to each reaction and incubated for 5 min to separate Tn5 from the DNA. Tagmentation purification with 1 × SPRI beads allowed for molecular barcode addition and PCR amplification to proceed to create NGS libraries. PCR reactions contained 1 unit of Kapa HiFi Polymerase (Kapa Biosystems), 1x HiFi buffer, 375 μM dNTPs, 0.5 μM of each primer, and purified tagmented sample. Cycling was initiated with gap-filling at 72°C for 3 min, followed by 98°C for 30 s. The PCR products were subjected to 10 cycles of denaturation at 63° C. for 30 s, annealing at 63° C. for 30 s, and extension at 72° C. for 3 min. Purification of the NGS libraries was performed using 1× SPRI beads.
[0128]
[0139] For patient samples, LR-PCR libraries were sequenced on a HiSeq 2500 (Illumina) in rapid run mode (paired-end, 150 cycles each). For cell line samples, LR-PCR libraries were sequenced on a NextSeq 550 (Illumina) to a minimum depth of 500 reads (single-end, 150 cycles).
[0129] Hybrid capture and sequencing:
[0140] Targeted NGS was performed as previously described [7, 8]. Briefly, DNA was isolated from patients' blood or saliva samples, quantified by a dye-based fluorescent assay, and then fragmented to 200–1000 bp by sonication. The fragmented DNA was converted into an NGS library by end-modification, A-tailing, and adapter ligation. Samples were then amplified by PCR with barcoded primers, multiplexed, and subjected to hybrid capture-based enrichment using 40-mer oligonucleotides (Integrated DNA Technologies) complementary to the common region between PMS2 and PMS2CL. NGS was performed on a HiSeq 2500 with an average sequencing depth of ∼500× for the entire panel (coverage at PMS2 was ∼1000×). All targeted nucleotides should be covered at a minimum depth of 20 reads.
[0130] Read alignment:
[0141] For hybrid capture data, paired-end NGS reads were first aligned to the hg19 human reference genome using BWA-MEM
[27] to tally reads originating from PMS2 and PMS2CL at the PMS2 locus of the reference genome. Alignments at exon 11 of PMS2 were filtered to include only reads overlapping at the site of known differences between the gene and pseudogene. Reads aligned to exons 12–15 of PMS2 and reads aligned to exons 3–6 of PMS2CL were partitioned into BAM files using samtool
[28] . The BAM files were converted into two unaligned FASTQ files (each number of parsed read pairs into one of two files) using Picard (Broad Institute). Each single-end FASTQ file was realigned separately against the hg19 genome to allow for ambiguous alignments, and reporting of several top alignments for each read. The resulting single-end alignments were used to generate paired-end alignments in the following manner: 1) both single-end reads had the same read name, 2) both single-end reads were mapped to a region spanning exons 12-15 of PMS2, 3) both single-end reads were aligned within 1000 bp of each other, and 4) a large number of putative pairs met the above criteria for a given read name, and the pair with the highest alignment score was selected. Reads unable to form suitable pairs as above were discarded. The resulting paired-end BAM file contained reads originating from both PMS2 and PMS2CL mapped against the PMS2 sequence.
[0131]
[0142] For the RT-PCR data (described below) and LR-PCR data, NGS reads were aligned to the hg19 genome sequence with the PMS2CL sequence removed, thereby tabulating gene and pseudogene reads in PMS2.
[0132] SNV and indel calling:
[0143] In the PMS2 region, where reads from PMS2 and PMS2CL were mapped (see above), SNVs and short indels were identified using GATK 4.0 HaplotypeCaller
[29] with sample ploidy option set to 4, max-reads-per-alignment-start option turned off, and min-pruning option set to 1. In the diploid PMS2 exon 11 region, SNVs and short indels were identified using GATK 1.6
[30] and FreeBayes
[31] . For diploid SNV calling in LR-PCR data, GATK 1.6 was used as well. In LR-PCR samples where we suspected allele dropout (see Discussion), AB was determined by visual inspection of the NGS data in the Integrative Genomics Viewer
[32] .
[0133] Calling CNV:
[0144] For short-read NGS of hybrid capture fragments, CNVs in exon 11 of PMS2 were determined by measuring the relative NGS read depth at the target location using a previously described algorithm [7]. To call CNVs in exons 12-15 of PMS2 from BAM files where reads originating from PMS2 and PMS2CL were located in the PMS2 sequence (see "Read alignment" above), two modifications to the CNV calling algorithm were made: 1) the predicted wild-type copy number was changed from 2 to 4 copies, and 2) P, a parameter that determines how likely the HMM will transition from the wild-type to the CNV state, was changed. CNV was set to 0.01, yielding high CNV sensitivity and specificity from empirical data.
[0134]
[0145] For CNV calls from LR-PCR data, read depths were counted in equal-sized bins (50 bp) flanking the amplicons. Bin counts for each sample were normalized by the median bin depth of the sample, and then each bin value was normalized by the bin median. The same bins were used for the corresponding regions of PMS2 and PMS2CL. The resulting binned and normalized data were searched for CNVs using a previously described algorithm [7]. Missing CNV calls were manually reviewed and status elucidated as positive or negative.
[0135] Simulation of CNV:
[0146] Single copy duplications and deletions were introduced by modifying the number of reads observed in one of the CNV-negative samples of a given batch of samples, as previously described
[33] . In exons 12–15 of PMS2, which had a baseline copy number of 4, single copy deletions and duplications were introduced by subsampling the reads by 75% or increasing the number of reads by 125%, respectively. Simulated CNVs were created for all possible consecutive combinations of exons in the four final exons of PMS2. For each CNV size and location, 2186 samples were simulated and tested by the CNV calling algorithm, and sensitivity was calculated as the percentage of correctly detected synthetic CNVs. Because pseudogene reads were filtered from the gene sequence, CNVs were simulated separately in exon 11 of PMS2, which had a baseline copy number of 2.
[0136] Simulation of tetraploid indels:
[0147] Indels in a tetraploid background (related to exons 12-15 of PMS2, where reads originating from genes and pseudogenes were remapped) were simulated to better test the sensitivity of indel calling using GATK4. Two diploid algorithms, at least one of which was previously determined by the Counsyl Reliant HCS panel to contain indels, were used. The diploid algorithm was merged to create a tetraploid alignment.If one of the samples has more reads than the other samples in the 100bp region that is located in the center of the indel, the reads were downsampled by binomial so that each merged diploid sample has approximately the same number of aligned reads.Then, GATK4 was used to call indels from these synthetic tetraploid alignments as described in the above section SNV and indel calling.
[0137] Variant Selection
[0148] For all variants in the five final exons of PMS2, variant interpretation was performed according to the American College of Medical Genetics and Genomics (ACMG) criteria, which use a five-tier classification category system (benign, likely benign, variant of uncertain pathogenic significance, likely pathogenic, pathogenic)
[34] . Classification was performed using available evidence in published literature and publicly available databases. Rules based on allele frequency were not used because identification of PMS2 variants in population databases may be inaccurate. Variant classification was reviewed and approved by board-certified laboratory directors.
[0138] MLPA:
[0149] MLPA was performed according to the manufacturer's protocols (MRC Holland, probemix P008-C1 PMS2 protocol published 12 / 11 / 17 and MLPA General Protocol published 3 / 23 / 18). In general, genomic DNA was covered with mineral oil to reduce evaporation during hybridization and ligation, then DNA was denatured at 98°C for 5 min and then held at 25°C. Hybridization reagent and probe mix were added to the samples and incubated at 95°C for 1 min and then at 60°C for 16-20 h. Probe pairs binding to closely spaced target DNA were ligated at 54°C for 15 min and then amplified by PCR for 35 cycles. Amplified probes were mixed with ROX ladder and formamide and then separated in a capillary electrophoresis instrument. PMS2 probe intensities were normalized against reference probe intensities, first within each sample and then between samples, by Coffalyser software (MRC Holland). The normalized probe intensity of each sample was compared to the average intensity of the reference samples, and Coffalyser made a CNV call in that region.
[0139] Reflex Rate Rating:
[0150] We used SNV-, indel-, and CNV-specific reflex rates from the LR-PCR and hybrid capture data and then extrapolated to large cohort sizes using Markov Chain Monte Carlo simulations with pymc
[35] to estimate reflex rates.
[0140] Identification of base analysis:
[0151] NGS reads from LR-PCR amplicons from PMS2 and PMS2CL were aligned to PMS2 and variants were called using the GATK Universal Genotyper. A site was considered reliable if the variant was homozygous to the reference allele in the PMS2-specific amplicon and homozygous to the alternative allele in the PMS2CL-specific amplicon (as aligned to PMS2) in 100% of the samples.
[0141] RNA testing: RNA extraction and reverse transcription:
[0152] RNA was extracted from 33 samples using the Agencourt RNAdvance Blood kit (Beckman Coulter) from 400 μL of whole blood according to the manufacturer's instructions. RNA was extracted from blood tubes within 7 days after blood collection. Extraction quality was assessed with the RNA 6000 Nano kit (Agilent). RNA was quantified with the Qubit HS RNA Assay kit (Thermo Fisher).
[0142]
[0153] RNA was reverse transcribed using Superscript II Reverse Transcriptase with oligo-dT and random hexamers as primers (kit from Thermo Fisher). Reactions were performed as follows: 0.1-1.0 μg total RNA, 1.25 μM of both random hexamer and oligo-dT primers, 0.8 mM dNTPs, and water in a final volume of 12 μL. Reactions were heated to 65°C for 5 min and then chilled on ice for 5 min. 1× First-Strand Buffer and 0.01 M DTT were added to each reaction and incubated at 42°C for 2 min. 10 U / μL Superscript II Reverse Transcriptase was added to each reaction and incubated at 42°C for 50 min and then heat inactivated at 72°C for 15 min. A positive control of pooled mRNA (Stratagene, catalog no. 750500-41) was used in each reverse transcription reaction.
[0143]
[0154] After reverse transcription, RNA was hydrolyzed with 2 μL of 1N NaOH and heated at 95°C for 5 min. Reactions were neutralized for downstream processing using 4 μL of 1M Tris-HCL (pH 7.5). cDNA was quantified using the Qubit ssDNA Assay kit (Thermo Fisher).
[0144] PCR:
[0155] For each sample, two reactions were set up: 1) forward primer PMS2_RNA_F and reverse primer RNA_Unv_R to amplify 1.5 kb of PMS2 from cDNA, and 2) forward primer PMS2CL_F and reverse primer RNA_Unv_R to amplify 1.5 kb of PMS2CL from cDNA (primer sequences in Table S3 in the Appendix). PCR reactions were performed in 1x LongAmp The PCR mixture contained Taq Reaction Buffer (NEB), 0.3 mM dNTPs, 1 μM each of forward and reverse primers, 20–70 ng of cDNA, 0.1 U / μL LongAmp Taq DNA polymerase (NEB), and water to 25 μL. Thermal cycling was as follows: 94°C for 5 min, 30 cycles of 94°C for 30 s, annealing at 52°C for PMS2 and 55°C for PMS2CL, 65°C for 2 min, followed by a final extension at 65°C for 10 min, then a 4°C hold. PCR products were purified with 1.2× SPRI beads. Amplicons were visualized on a 2% agarose gel or with a DNA7500 kit (Agilent).
[0145] Sequencing:
[0156] 50–100 ng of each amplicon was sheared in a 50 μL volume using a Bioruptor (Diagenode) with 12 cycles of 30 s on and 90 s off. Fragmentation was visualized with a High Sensitivity DNA kit (Agilent). All sheared material was used as input for library preparation. The KAPA Hyper Prep kit (Kapa Biosystems) was used for library preparation and the manufacturer's instructions were followed. Adaptors were diluted to 15 μM for PMS2 and 3 μM for PMS2CL. Enrichment PCR was performed for 9 cycles. Absorbance measurements (Tecan Samples were quantified using a 500 mM M200 kit, normalized to 10 nM, and standardized to one reaction. Final libraries were quantified by qPCR using the KAPA Library Quantification Kit (Kapa Biosystems) and double indexed. The single reads containing the DNA were sequenced for 75 cycles using a NextSeq 550 System (Illumina).
[0146] alignment:
[0157] The base call files were converted to FASTQ files using bcl2fastq (Illumina), and the FASTQ files were aligned using STAR
[36] .
[0147] Analysis metrics:
[0158] Metrics were defined as follows: Sensitivity = TP / (TP+FN); Specificity = TN / (TN+FP). CIs were calculated by the method of Clopper and Pearson
[37] . For SNVs and indels, true negatives were defined as concordant negative results observed at sites found to be polymorphic in the cohort used (positions where we observed a non-canonical base in at least one sample).
[0148] result Zero nucleotides can reliably identify exons 12-15 of PMS2 as PMS2CL:
[0159] NGS of short DNA fragments will be able to identify PMS2-specific variants in the five final exons only if the fragments themselves can be unambiguously aligned to the gene or pseudogene. To overcome the interference of pseudogenes, unique mapping will rely on bases that differ between PMS2 and PMS2CL. In the hg19 reference genome, these discriminating bases are rare (Figure 1D, left bar): sequence identity in each of the five final exons of PMS2 (filled in with 20 nt of intronic sequence) is over 97%, with differences including only 26, 0, 1, 1, and 0 bases in exons 11 to 15, respectively. Furthermore, previous reports have shown that natural mutations may suppress the reliability of these discriminating bases represented in the reference genome [17, 18].
[0149]
[0160] To test the reliability of the reference genome, we assembled a set of natural variants in exons 11-15 of PMS2 and the corresponding regions in PMS2CL. NGS was performed on gene- and pseudogene-specific LR-PCR amplicons for 707 patient samples of the used cohort (Table 1) with various self-reported ethnic affiliations (Appendix Table S4). Seven of the 26 predicted positions in exon 11 of PMS2 were found to have distinct alleles in the gene and pseudogene, making them reliable discriminators. In contrast, for 19 positions in exon 11 and two positions in exons 12-15, alleles ostensibly specific to PMS2 from hg19 were observed at least once in the PMS2CL LR-PCR data and vice versa (see Appendix Table S4 for allele frequencies). Thus, after taking into account natural variation in the gene and pseudogenes, there are no reliable discriminatory bases in exons 12–15 of PMS2 (i.e., 100% sequence identity) and seven discriminatory bases in exon 11 (Figure 1D, dark bars). Collectively, these data suggest that variant identification by short-read NGS alone is sufficient for exon 11, but that a different approach is required for exons 12–15.
[0150] [Table 1]
[0151] Reflex workflow for disambiguating variants discovered in short-read NGS:
[0161] Using short-read NGS as the basis, we assessed the validity of a workflow for the 3' exon of PMS2, performing reflex testing, including orthogonal assays to clarify whether variants are of genic or pseudogenic origin, only if clinically indicated (Figure 2A). In the short-read NGS stage of the study, the molecular approach is concordant across the five final exons of PMS2. By designing capture probes that specifically avoid positions shown to be altered between PMS2 and PMS2CL in LR-PCR data from patient samples, we capture DNA fragments in a manner that is ambiguous as to whether they are of genic or pseudogenic origin (Figure 2B, purple box).
[0152]
[0162] The workflow uses different bioinformatics strategies for exon 11 and exons 12-15 of PMS2 (Figure 2B, blue boxes). For exon 11, PMS2-specific variants are identified by adjusting the read-alignment software to partition reads relative to PMS2 or PMS2CL based on gene and pseudogene discriminant bases. In contrast, for exons 12-15 of PMS2, reads are aligned with permissive settings such that each read aligns to its best gene position and its best pseudogene position (see Methods). In a typical sample with two copies each for PMS2 and PMS2CL, this approach effectively results in a read depth at each position corresponding to four copies. To identify SNVs, indels, and CNVs, the variant calling software is adjusted to expect a baseline ploidy of two in exon 11 and four in exons 12-15 (Figure 2B, blue and green boxes).
[0153]
[0163] Disambiguation by reflex testing is only required for a subset of variants based on their type and clinical interpretation (Figure 2B, orange box). Thus, variant interpretation is performed prior to reflex testing. Benign variants are not reflex tested or reported to patients. Samples with CNVs in any of the five final exons of PMS2 classified as pathogenic, likely pathogenic, or variants of unknown pathogenic significance (VUS) undergo reflex testing for disambiguation. Samples containing non-benign SNVs or indels in exons 12-15 are reflex tested for disambiguation, whereas samples with such variants in exon 11 are only reported without reflex by unique read mapping within that exon. Disambiguation testing for SNVs, indels, and CNVs can be performed by LR-PCR followed by sequencing to determine whether the variant is from PMS2 or PMS2CL; MLPA can assist in resolving CNVs
[20] .
[0154]
[0164] By implementing the proposed workflow, the cancer risk associated with the five final exons of PMS2 is elucidated for the majority of samples using only short-read NGS. For each of the 707 patient samples that underwent LR-PCR (Table 1), variant classification was performed on the results, and it was found that almost 93% did not require reflex testing. The remaining approximately 7% required subsequent testing to obtain a confident PMS2 screening result (Figure 2A). The SNV and indel specific components of this reflex rate were 41 / 707 (5.8%), while the reflex rates due to CNV calls and no calls were 2 / 707 (0.3%) and 1 / 144 (0.7%), respectively. Using simulations (see Methods), we estimated the reflex rate for a large cohort of 13,000 patients to be 7.7% (95% CI: 5.4-10.7%). The 0.7% contribution to the reflex rate from samples without CNV calls is predicted to be an upper-bound estimate, since the standard practice of retesting such samples at least once on short-read NGS yields confident negative calls (data not shown), thereby avoiding reflex testing. Thus, the overall reflex rate of the proposed workflow (see Figure 6) is expected to be less than 8%.
[0155] Short-read NGS accurately identified samples requiring reflex testing for SNVs and indels:
[0165] The reflex workflow described herein is clinically feasible only if short-read NGS testing (Figure 2) has high analytical sensitivity and specificity to (1) identify variants in exon 11 of PMS2 and (2) signal samples requiring reflex testing for variants in exons 12-15 with ambiguity of PMS2 / PMS2CL origin. To evaluate the accuracy of short-read NGS testing for SNVs and indels, we compared the results to those observed with LR-PCR for 144 patient samples and 155 cell lines (Figure 3). An irregular confusion matrix was required by measuring genotype concordance in exons 12-15 because short-read NGS genotypes are reported as tetraploid (see Methods), whereas LR-PCR returns diploid genotype calls for both genes and pseudogenes (Figure 3A highlights some examples). The matrix includes "tolerated dosage errors" where the presence of alternative alleles is adequately detected but the number of alternative alleles is not consistent; such errors are considered tolerable because the presence of alternative alleles in short-read NGS is sufficient to trigger reflex testing and be corrected. When compared at 1,678 sites using LR-PCR as the true set, the short-read NGS test had 100% analytical sensitivity and 100% analytical specificity in exon 11 (Figure 3B) and 99.9% analytical sensitivity and 100% analytical specificity in exons 12 to 15 (Figure 3C).
[0156]
[0166] The paucity of indel calls in the patient cohorts and cell lines used (17 in total), coupled with the rare use of variant calling software in tetraploid-background mode for clinical genomic applications, motivated a deeper investigation of indel calling efficiency in exons 12–15 of PMS2. Expected NGS data were simulated for samples with a tetraploid genomic background populated with indels of different allele dosages (1, 2, 3, or 4 copies). To construct such samples, To do this, we merged diploid NGS data from two samples (at least one containing an indel) in the region of the HCS test used outside of PMS2 (Figure 4A, see Methods). The genotypes of each of the two samples yield the predicted genotype of the merged sample, e.g., combining a homozygous alternative sample (two indel alleles) with a heterozygous sample (one indel allele) would give a predicted indel dosage of 3. Figure 4B shows a sensitivity of 99.6% for indels in a simulated tetraploid background, suggesting that sensitivity is relatively high in exons 12-15 of PMS2, where the read alignment and variant calling strategy used results in a tetraploid background. Because the empirical data in Figure 3C demonstrate 100% specificity for indels in exons 12-15, specificity was not further evaluated in the simulations.
[0157]
[0167] Taken together, comparison of SNV and indel calling between LR-PCR and short-read NGS suggests that the pre-reflex step of the proposed workflow described herein achieves sufficient analytical sensitivity and specificity to be considered for clinical use.
[0158] Accurate detection of short-read NGS for samples requiring CNV reflex testing
[0168] To evaluate the sensitivity and specificity of short-read NGS for CNVs in the five final exons of PMS2, patient samples, cell lines, known positives, and samples with simulated positives were tested. Similar to SNVs and indels, the CNV detection algorithm described above was adapted to use copy number baselines of two for exon 11 and four in exons 12-15 of PMS2 (Figure 2B, blue boxes; see Methods). Three known positive samples with CNVs in the five final exons were correctly identified as having CNVs encompassing the predicted exons (Figure 5A). We further observed a deletion of exons 13-14 in four of the cell lines and one of the clinical samples; in the clinical sample, short-read NGS identified a signal drop from the tetraploid background (Figure 5B), MLPA confirmed the presence of a similar deletion (Figure 5C), and NGS on LR-PCR amplicons revealed that the deletion was present in PMS2CL rather than PMS2 (Figure 5D). Interestingly, although only one of the two copies of this region is deleted in PMS2CL, the LR-PCR profile shows a 75% signal drop in the deleted region. During LR-PCR, this is speculated to result from preferential amplification of the allele carrying the shorter deletion. Thus, the LR-PCR data was unique in providing disambiguation, whereas the short-read NGS and MLPA data had more easily interpretable copy number values.
[0159]
[0169] Due to the absence of a large set of CNV-positive samples, a complete and straightforward characterization of the sensitivity of PMS2 CNV calling on short-read NGS would require blinded testing of thousands of samples. Instead, sequencing data from a large number of CNV-negative patients was used as a substrate in simulations introducing CNVs of a given length and position (see Methods). By running the CNV detection algorithm described above on 2186 simulated samples, we measured the analytical sensitivity for CNVs ranging from 1 to 5 exon lengths (Table 2; simulated data for cell line samples in Appendix Table S6). Overall, the sensitivity for multiple exon deletions was greater than 99.2%, and the sensitivity for single exon deletions was approximately 89%. By weighting the simulated sensitivities by the observed frequency distribution of CNV lengths in the five final exons of PMS2 [21, 23, 24], the aggregate CNV sensitivity in this complex genomic region is estimated to be 96.7%.
[0160] [Table 2]
[0161]
[0170] High sensitivity for CNV does not come at the expense of low specificity. This triggers the measurement of the CNV false positive rate for the large cohort used. In the 302 hybrid capture cohort of 302 samples, there was one no call, which was treated as a false positive. Thus, the sample-level specificity is 99.7% (95% CI: 98.2-100%).
[0162]
[0171] Based on these analyses, it was concluded that short-read NGS (optimized in the described workflow) could achieve >96% sensitivity and >99% specificity for detecting samples containing CNVs in the five terminal exons of PMS2. Gene and pseudogene specific variant information for common cell lines:
[0163]
[0172] Reference cell lines with known genotypes facilitate the development and evaluation of novel molecular diagnostic methods, but samples with high-quality genotypes in the PMS2 region are generally unavailable due to the complex nature of the region. In the course of developing and testing the workflow characterized above, we performed NGS on both hybrid capture fragments and LR-PCR amplicons in cell lines assembled from whole genome sequencing (Illumina Polaris 1 Diversity Panel) or Genome in the Bottle (GIAB) Consortium with high-quality genome sequences at a depth of approximately 30× [38,39]. Importantly, Figure 7 shows that the observed gene-specific genotypes differed from the Polaris and GIAB data (including phase data for GIAB samples; Figure 7C). In principle, such differences could arise in part due to errors in either dataset, e.g., due to biological contamination, non-specific amplification, non-specific sequence alignment, or technical processing errors by the selected genotyping software. Although the agreement between the orthogonal hybrid capture and LR-PCR assays suggests that the genotypes reported in this invention are correct, as a third orthogonal method, genotyping of PMS2 and PMS2CL was performed from RNA extracted from 33 of the LR-PCR samples (see Methods). The RNA-derived genotypes were consistent with the LR-PCR data (Figure 8), strongly suggesting that we have defined the correct gene- and pseudogene-specific genotypes. Gene- and pseudogene-specific variant information will be shared to aid scientific research and clinical development of PMS2 and its role in Lynch syndrome. For patient samples, variant frequencies will be given (Appendix Table S4) to share valuable data, while keeping in mind patient consent and PHI compliance. For cell lines, variant frequencies and BAM and VCF files for LR-PCR amplicons spanning the five final exons of PMS2 and PMS2CL will be shared. (Appendix Table S5 and ENA accession no. PRJEB27948).
[0164] Exemplary embodiments
[0173] The following embodiments are illustrative and are not intended to limit the present invention.
[0165]
[0174] Embodiment 1. A method for detecting genetic variations in a genome of a subject, the genome comprising a first region and a second region of interest that are highly homologous, the method comprising: (a) obtaining sequence reads by paired-end sequencing from multiple sites of a target in a first region and a second region of the target, the sequence reads including a first read and a second read obtained at each site of the target; (b) aligning the sequence reads to a reference genome, wherein the first read and the second read are aligned separately to the reference genome, and the aligner generates multiple possible alignments for each of the first read and the second read; (c) identifying a first read and a second read that align to a first region of the object; (d) pairing a first read and a second read from the reads identified in step (c), thereby generating a top pair alignment; (e) detecting genetic variations in the top pair alignment generated in step (d); A method comprising:
[0166]
[0175] Embodiment 2. The method of embodiment 1, comprising, prior to step (b), a step of aligning the first read and the second read to a reference genome, wherein for each pair of the first read and the second read, the aligner generates a best possible paired-end alignment for the first or second region of the object, and only the paired-end reads associated with the top alignment scores for the first or second region of the object are aligned separately in step (b).
[0167]
[0176] Embodiment 3. The method of embodiment 1, wherein the sequence reads are obtained by direct targeted sequencing (DTS) of multiple sites of interest, a first read comprising a genomic sequence read and a second read comprising a probe sequence read associated with the sites of interest.
[0168]
[0177] Embodiment 4. The method of embodiment 1, wherein in step (b), the sequence reads are aligned using the Burrows-Wheeler Aligner (BWA) algorithm.
[0169]
[0178] Embodiment 5. The method of embodiment 1, wherein in step (b), the aligner only issues alignments that meet a minimum alignment score for the first and second regions of the object.
[0170]
[0179] Embodiment 6. The method of embodiment 1, wherein the first read and the second read are paired in step (d) only if the alignment of the first read and the second read to the first region of interest is within a certain number of bases of each other.
[0171]
[0180] Embodiment 7. The alignment of the first read and the second read to the first region of the object is about 100 bp, about 200 bp, about 200 bp, about 300 bp, about 400 bp, about 500 bp, about 600 bp, about 700 bp, about 800 bp, about 900 bp, about 1000 bp, about 1100 bp, about 1200 bp, about 1300 bp, about 1400 bp, about 1500 bp, about 1600 bp, about 1700 bp, about 1800 bp, about 1900 bp, about 2100 bp, about 2200 bp, about 2300 bp, about 2400 bp, about 2500 bp, about 2600 bp, about 2700 bp, 2. The method of embodiment 1, wherein the first read and the second read are paired in step (d) only if they are within 500 bp, or more than 1500 bp.
[0172]
[0181] Embodiment 8. The method of embodiment 1, comprising in step (d) generating a number of pairwise alignments, calculating an alignment score for each of the number of pairwise alignments, and identifying a top pairwise alignment having the highest alignment score.
[0173]
[0182] Embodiment 9. The method of embodiment 1, wherein the top pair alignment in step (d) is selected as the one having the smallest template length.
[0174]
[0183] Embodiment 10. The method of embodiment 1, wherein the genetic variation comprises a SNP, an indel, an inversion, and / or a CNV.
[0175]
[0184] Embodiment 11. The method of embodiment 1, wherein the detecting in step (e) comprises calling SNPs, indels, inversions, and / or CNVs.
[0176]
[0185] Embodiment 12. The method of embodiment 1, wherein the detecting in step (e) comprises using a Hidden Markov Model (HMM) caller to determine copy number.
[0177]
[0186] Embodiment 13. The method of embodiment 1, wherein the detecting in step (e) is based on a predicted ploidy of 2.
[0178]
[0187] Embodiment 14. The method of embodiment 1, wherein the detecting in step (e) is based on a predicted ploidy of 4.
[0179]
[0188] Embodiment 15. The method of embodiment 1, wherein if a genetic variation is detected in step (e), a portion of the subject's genome is amplified by long-range PCR and assayed by multiplex ligation-dependent probe amplification (MLPA).
[0180]
[0189] Embodiment 16. The method according to embodiment 1, wherein if a genetic variation is detected in step (e), a portion of the first region of interest is amplified by long-range PCR, and the product or a portion thereof is sequenced by Sanger sequencing or NGS.
[0181]
[0190] Embodiment 17. The method of embodiment 1, wherein if a genetic variation is detected in step (e), the subject's genomic DNA is assayed by multiplex ligation-dependent probe amplification (MLPA).
[0182]
[0191] Embodiment 18. The method of embodiment 1, wherein the sequence reads are 30-50 bp or 100-200 bp in length.
[0183]
[0192] Embodiment 19. The method of embodiment 1, wherein the highly homologous first and second regions of interest are at least 80%, at least 81%, at least 82%, at least 83%, at least 84%, at least 85%, at least 86%, at least 87%, at least 88%, at least 89%, at least 90%, at least 91%, at least 92%, at least 93%, at least 94%, at least 95%, at least 96%, at least 97%, at least 98%, at least 99%, or greater than 99% identical.
[0184]
[0193] Embodiment 20. The method of embodiment 1, wherein the sequence reads are obtained from one or more exons within the first region and / or the second region of interest.
[0185]
[0194] Embodiment 21. The method of embodiment 1, wherein the sequence reads are obtained from one or more introns within the first region and / or the second region of interest.
[0186]
[0195] Embodiment 22. The method of embodiment 1, wherein the sequence reads are obtained from one or more exons and introns within the first region and / or the second region of interest.
[0187]
[0196] Embodiment 23. The method of embodiment 1, wherein sequence reads are obtained from one or more exons and introns within the first region and / or second region of interest, and the introns are located adjacent to the exons.
[0188]
[0197] Embodiment 24. The method of embodiment 1, wherein the sequence reads are obtained from one or more clinically addressable regions associated with the first region and / or the second region of interest.
[0189]
[0198] Embodiment 25 The method of embodiment 1, wherein the first region of interest comprises a gene and the second region of interest comprises a pseudogene.
[0190]
[0199] Embodiment 26 The method of embodiment 1, wherein the first region of interest comprises a pseudogene and the second region of interest comprises a gene.
[0191]
[0200] Embodiment 27 The method of embodiment 1, wherein the first region of interest comprises two alleles.
[0192]
[0201] Embodiment 28 The method of embodiment 1, wherein the second region of interest comprises two alleles.
[0193]
[0202] Embodiment 29. The method of any one of embodiments 25 to 28, wherein the gene is PMS2.
[0194]
[0203] Embodiment 30. The method of any one of embodiments 25 to 28, wherein the pseudogene is PMS2CL.
[0195]
[0204] Embodiment 31. The method of embodiment 1, wherein multiple sites of interest are present within an exon of PMS2 and an exon of another portion of the subject's genome.
[0196]
[0205] Embodiment 32. The method of embodiment 1, wherein the multiple sites of interest are present within an exon of PMS2 and an exon of PMS2CL.
[0197]
[0206] Embodiment 33. The method of embodiment 1, wherein multiple sites of interest are present within exons 11, 12, 13, 14, and / or 15 of PMS2 and exons 2, 3, 4, 5, and / or 6 of PMS2CL.
[0198]
[0207] Embodiment 34 The method of embodiment 1, wherein the subject is a human and the sequence reads are aligned to a human reference genome.
[0199]
[0208] Embodiment 35. The method of embodiment 1, implemented by a computer.
[0200]
[0209] Embodiment 36. The method of embodiment 1, wherein the reference genome does not contain a masked or modified portion of the first homologous region or the second homologous region of interest.
[0201]
[0210] Embodiment 37. A non-transitory computer-readable storage medium containing computer-executable instructions for implementing embodiment 1.
[0202]
[0211] 38. (a) one or more processors; (b) memory, and (c) one or more programs A system comprising: one or more programs stored in a memory and configured to be executed by one or more processors, the one or more programs including instructions for carrying out embodiment 1.
[0203]
[0212] References 1. Nagy R, Sweet K, Eng C. Highly penetrant hereditary cancer syndromes. Oncogene. 2004;23: 6445-6470. 2. Lu KH, Wood ME, Daniels M, Burke C, Ford J, Kauff ND, et al. American Society of Clinical Oncology Expert Statement: collection and use of a cancer family history for oncology providers. J Clin Oncol. 2014;32: 833-840. 3. Mucci LA, Hjelmborg JB, Harris JR, Czene K, Havelick DJ, Scheike T, et al. Familial Risk and Heritability of Cancer Among Twins in Nordic Countries. JAMA. 2016;315: 68-76. 4. Foulkes WD. Inherited Susceptibility to Common Cancers. N Engl J Med. 2008;359: 2143-2153. 5. Garber JE, Offit K. Hereditary cancer predisposition syndromes. J Clin Oncol. 2005;23: 276-292. 6. Vogelstein B, Papadopoulos N, Velculescu VE, Zhou S, Diaz LA, Kinzler KW. Cancer Genome Landscapes. Science. 2013;339: 1546-1558. 7. Vysotskaia VS, Hogan GJ, Gould GM, Wang X, Robertson AD, Haas KR, et al. Development and validation of a 36-gene sequencing assay for hereditary cancer risk assessment. PeerJ. 2017;5: e3046. 8. Kang HP, Maguire JR, Chu CS, Haque IS, Lai H, Mar-Heyming R, et al. Design and validation of a next generation sequencing assay for hereditary BRCA1 and BRCA2 mutation testing. PeerJ. 2016;4: e2162. 9. Bunnell AE, Garby CA, Pearson EJ, Walker SA, Panos LE, Blum JL. The Clinical Utility of Next Generation Sequencing Results in a Community-Based Hereditary Cancer Risk Program. J Genet Couns. 2017;26: 105-112. 10. Desmond A, Kurian AW, Gabree M, Mills MA, Anderson MJ, Kobayashi Y, et al. Clinical Actionability of Multigene Panel Testing for Hereditary Breast and Ovarian Cancer Risk Assessment. JAMA Oncol. 2015;1: 943-951. 11. Lynch HT, Smyrk T, Lynch J, Fitzgibbons R Jr, Lanspa S, McGinn T. Update on the differential diagnosis, surveillance and management of hereditary non-polyposis colorectal cancer. Eur J Cancer. 1995;31A: 1039-1046. 12. Blount J, Prakash A. The changing landscape of Lynch syndrome due to PMS2 mutations. Clin Genet. 2018;94: 61-69. 13. Sijmons RH, Hofstra RMW. Review: Clinical aspects of hereditary DNA Mismatch repair gene mutations. DNA Repair . 2016;38: 155-162. 14. Tiwari AK, Roy HK, Lynch HT. Lynch syndrome in the 21st century: clinical perspectives. QJM. 2016;109: 151-158. 15. Lynch HT, Fusaro RM, Lynch JF. Cancer Genetics in the New Era of Molecular Biology. Ann N Y Acad Sci. 1997;833: 1-28. 16. De Vos M, Hayward BE, Picton S, Sheridan E, Bonthron DT. Novel PMS2 pseudogenes can conceal recessive mutations causing a distinctive childhood cancer syndrome. Am J Hum Genet. 2004;74: 954-964. 17. Hayward BE, De Vos M, Valleley EMA, Charlton RS, Taylor GR, Sheridan E, et al. Extensive gene conversion at the PMS2 DNA mismatch repair locus. Hum Mutat. 2007;28: 424-430. 18. van der Klift HM, Tops CMJ, Bik EC, Boogaard MW, Borgstein A-M, Hansson KBM, et al. Quantification of sequence exchange events between PMS2 and PMS2CL provides a basis for improved mutation scanning of Lynch syndrome patients. Hum Mutat. 2010;31: 578-587. 19. Vaughn CP, Robles J, Swensen JJ, Miller CE, Lyon E, Mao R, et al. Clinical analysis of PMS2: mutation detection and avoidance of pseudogenes. Hum Mutat. 2010;31: 588-593. 20. Vaughn CP, Hart KJ, Samowitz WS, Swensen JJ. Avoidance of pseudogene interference in the detection of 3' deletions in PMS2. Hum Mutat. 2011;32: 1063-1071. 21. van der Klift HM, Mensenkamp AR, Drost M, Bik EC, Vos YJ, Gille HJJP, et al. Comprehensive Mutation Analysis of PMS2 in a Large Cohort of Probands Suspected of Lynch Syndrome or Constitutional Mismatch Repair Deficiency Syndrome. Hum Mutat. 2016;37: 1162-1179. 22. Li J, Dai H, Feng Y, Tang J, Chen S, Tian X, et al. A Comprehensive Strategy for Accurate Mutation Detection of the Highly Homologous PMS2. J Mol Diagn. 2015;17: 545-553. 23. Vaughn CP, Baker CL, Samowitz WS, Swensen JJ. The frequency of previously undetectable deletions involving 3’ Exons of the PMS2 gene. Genes Chromosomes Cancer. 2013;52: 107-112. 24. Espenschied CR, LaDuca H, Li S, McFarland R, Gau C-L, Hampel H. Multigene Panel Testing Provides a New Perspective on Lynch Syndrome. J Clin Oncol. 2017;35: 2568-2575. 25. Etzler J, Peyrl A, Zatkova A, Schildhaus H-U, Ficek A, Merkelbach-Bruse S, et al. RNA-based mutation analysis identifies an unusual MSH6 splicing defect and circumvents PMS2 pseudogene interference. Hum Mutat. 2008;29: 299-305. 26. Herman DS, Smith C, Liu C, Vaughn CP, Palaniappan S, Pritchard CC, et al. Efficient Detection of Copy Number Mutations in PMS2 Exons with a Close Homolog. J Mol Diagn. 2018;20: 512-521. 27. Li H. Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM [Internet]. 2013. Available: arxiv.org / abs / 1303.3997 28. Li H, Handsaker B, Wysoker A, Fennell T, Ruan J, Homer N, et al. The Sequence Alignment / Map format and SAMtools. Bioinformatics. 2009;25: 2078-2079. 29. Poplin R, Ruano-Rubio V, DePristo MA, Fennell TJ, Carneiro MO, Van der Auwera GA, et al. Scaling accurate genetic variant discovery to tens of thousands of samples [Internet]. 2017. doi:10.1101 / 201178 30. Garrison E, Marth G. Haplotype-based variant detection from short-read sequencing [Internet]. arXiv [q-bio.GN]. 2012. Available: arxiv.org / abs / 1207.3907 31. McKenna A, Hanna M, Banks E, Sivachenko A, Cibulskis K, Kernytsky A, et al. The Genome Analysis Toolkit: a MapReduce framework for analyzing next-generation DNA sequencing data. Genome Res. 2010;20: 1297-1303. 32. Home | Integrative Genomics Viewer [Internet]. [cited 7 Sep 2018]. Available: www.broadinstitute.org / igv 33. Hogan GJ, Vysotskaia VS, Beauchamp KA, Seisenberger S, Grauman PV, Haas KR, et al. Validation of an Expanded Carrier Screen that Optimizes Sensitivity via Full-Exon Sequencing and Panel-wide Copy Number Variant Identification. Clin Chem. 2018;64: 1063-1073. 34. Richards S, Aziz N, Bale S, Bick D, Das S, Gastier-Foster J, et al. Standards and guidelines for the interpretation of sequence variants: a joint consensus recommendation of the American College of Medical Genetics and Genomics and the Association for Molecular Pathology. Genet Med. 2015;17: 405-424. 35. Salvatier J, Wiecki TV, Fonnesbeck C. Probabilistic programming in Python using PyMC3. PeerJ Comput Sci. PeerJ Inc.; 2016;2: e55. 36. Dobin A, Davis CA, Schlesinger F, Drenkow J, Zaleski C, Jha S, et al. STAR: ultrafast universal RNA-seq aligner. Bioinformatics. 2013;29: 15-21. 37. Clopper CJ, Pearson ES. The Use of Confidence or Fiducial Limits Illustrated in the Case of the Binomial. Biometrika. 1934;26: 404. 38. Zook JM, Catoe D, McDaniel J, Vang L, Spies N, Sidow A, et al. Extensive sequencing of seven human genomes to characterize benchmark reference materials. Sci Data. 2016;3: 160025. 39. Zook JM, Chapman B, Wang J, Mittelman D, Hofmann O, Hide W, et al. Integrating human sequence data sets provides a resource of benchmark SNP and indel genotype calls. Nat Biotechnol. 2014;32: 246-251.
[0213] It is understood that the examples and embodiments described herein are for illustrative purposes only, and that various modifications or changes in light thereof will be suggested to those skilled in the art and are to be included within the spirit and scope of this application and the scope of the appended claims. All publications, patents, and patent applications cited herein are incorporated herein in their entirety by reference for all purposes.
Claims
1. 1. A method for detecting genetic variations in a genome of a subject, the genome comprising first and second highly homologous regions of interest, the method comprising: (a) obtaining sequence reads by paired-end sequencing from multiple sites of a target in a first region and a second region of the target, the sequence reads including a first read and a second read obtained at each site of the target; (b) aligning the sequence reads to a reference genome, wherein the first read and the second read are aligned separately to the reference genome, and the aligner generates multiple possible alignments for each of the first read and the second read; (c) identifying a first read and a second read that align to a first region of the object; (d) pairing a first read and a second read from the reads identified in step (c), thereby generating a top pair alignment; (e) detecting genetic variations in the top pair alignments generated in step (d); The method comprising:
2. 2. The method of claim 1, further comprising, prior to step (b), a step of aligning the first read and the second read to a reference genome, wherein for each pair of the first read and the second read, the aligner generates a best possible paired-end alignment for the first or second region of the object, and only paired-end reads associated with the top alignment scores for the first or second region of the object are aligned separately in step (b).
3. 2. The method of claim 1, wherein the sequence reads are obtained by direct targeted sequencing (DTS) of multiple sites of interest, a first read comprising a genomic sequence read and a second read comprising a probe sequence read associated with the sites of interest.
4. 2. The method of claim 1, wherein in step (b), the sequence reads are aligned using the Burrows-Wheeler Aligner (BWA) algorithm.
5. The method of claim 1 , wherein in step (b), the aligner generates only alignments that satisfy a minimum alignment score for the first and second regions of the object.
6. 2. The method of claim 1, wherein the first read and the second read are paired in step (d) only if alignments of the first read and the second read to the first region of interest are within a certain number of bases of each other.
7. 2. The method of claim 1, wherein the first read and the second read are paired in step (d) only if the alignment of the first read and the second read to the first region of interest is within about 100 bp, about 200 bp, about 200 bp, about 300 bp, about 400 bp, about 500 bp, about 600 bp, about 700 bp, about 800 bp, about 900 bp, about 1000 bp, about 1100 bp, about 1200 bp, about 1300 bp, about 1400 bp, about 1500 bp, or more than 1500 bp.
8. 2. The method of claim 1, further comprising, in step (d), generating a number of pairwise alignments, calculating an alignment score for each of the number of pairwise alignments, and identifying a top pairwise alignment as having the highest alignment score.
9. The method of claim 1, wherein the top pair alignment in step (d) is selected as the one having the smallest template length.
10. The method of claim 1 , wherein the genetic variation comprises a SNP, an indel, an inversion, and / or a CNV.
11. 2. The method of claim 1, wherein the detecting in step (e) comprises calling SNPs, indels, inversions, and / or CNVs.
12. 10. The method of claim 1, wherein the detecting in step (e) comprises using a Hidden Markov Model (HMM) caller to determine copy number.
13. 2. The method of claim 1 , wherein the detecting in step (e) is based on a predicted ploidy of 2.
14. 2. The method of claim 1, wherein the detecting in step (e) is based on a predicted ploidy of 4.
15. 2. The method of claim 1, wherein if a genetic variation is detected in step (e), a portion of the subject's genome is amplified by long-range PCR and assayed by multiplex ligation-dependent probe amplification (MLPA).
16. 2. The method of claim 1, wherein if a genetic variation is detected in step (e), a portion of the first region of interest is amplified by long-range PCR and the product or a portion thereof is sequenced by Sanger sequencing or NGS.
17. 2. The method of claim 1, wherein if a genetic variation is detected in step (e), the subject's genomic DNA is assayed by multiplex ligation-dependent probe amplification (MLPA).
18. 2. The method of claim 1, wherein the sequence reads are 30-50 bp or 100-200 bp in length.
19. 2. The method of claim 1, wherein the highly homologous first and second regions of the target are at least 80%, at least 81%, at least 82%, at least 83%, at least 84%, at least 85%, at least 86%, at least 87%, at least 88%, at least 89%, at least 90%, at least 91%, at least 92%, at least 93%, at least 94%, at least 95%, at least 96%, at least 97%, at least 98%, at least 99%, or greater than 99% identical.
20. 2. The method of claim 1, wherein the sequence reads are obtained from one or more exons within the first region and / or the second region of interest.
21. The sequence reads identify one or more regions within the first and / or second region of interest. The method of claim 1, wherein the nucleic acid sequence is obtained from a number of introns.
22. 2. The method of claim 1, wherein the sequence reads are obtained from one or more exons and introns within the first region and / or the second region of interest.
23. 2. The method of claim 1, wherein sequence reads are obtained from one or more exons and introns within the first and / or second regions of interest, and the introns are located adjacent to the exons.
24. 2. The method of claim 1, wherein the sequence reads are obtained from one or more clinically addressable regions associated with the first region and / or the second region of interest.
25. 2. The method of claim 1, wherein the first region of interest comprises a gene and the second region of interest comprises a pseudogene.
26. 2. The method of claim 1, wherein the first region of interest comprises a pseudogene and the second region of interest comprises a gene.
27. The method of claim 1 , wherein the first region of interest comprises two alleles.
28. The method of claim 1 , wherein the second region of interest comprises two alleles.
29. The method according to any one of claims 25 to 28, wherein the gene is PMS2.
30. The method according to any one of claims 25 to 28, wherein the pseudogene is PMS2CL.
31. 2. The method of claim 1, wherein multiple sites of interest are present within an exon of PMS2 and an exon in another portion of the subject's genome.
32. 2. The method of claim 1, wherein multiple sites of interest are present within an exon of PMS2 and an exon of PMS2CL.
33. 2. The method of claim 1, wherein multiple sites of interest are present within exons 11, 12, 13, 14, and / or 15 of PMS2 and exons 2, 3, 4, 5, and / or 6 of PMS2CL.
34. 2. The method of claim 1, wherein the subject is a human and the sequence reads are aligned to a human reference genome.
35. The method of claim 1 , wherein the method is computer-implemented.
36. The method of claim 1 , wherein the reference genome does not contain masked or modified portions of the first homologous region or the second homologous region of interest.
37. 13. A non-transitory computer-readable storage medium comprising computer-executable instructions for carrying out claim 1.
38. (a) one or more processors; (b) memory; and (c) one or more programs 11. A system comprising: one or more programs stored in memory and configured to be executed by one or more processors, the one or more programs comprising instructions for carrying out claim 1.
Citation Information
Patent Citations
Nucleic acid sequence assembly
JP2018512092A
Methods, systems and processes of identifying genetic variation in highly similar genes
US20160300014A1