Unlocking Genetic Profiles with a Programmable DNA‐Powered Decoding Circuit
Introduction Genetic variations contribute dramatically to the inter‐individual phenotypic variability, ranging from eye color to disease susceptibility.  Over the past decade, tremendous advances in genomic technologies have revolutionized the ability of genome‐wide association studies (GWAS) to pinpoint human genetic variations implicated in complex diseases, such as cancer and type 2 diabetes (T2D).  Discoveries of numerous risk variants not only shed light on the pathogenesis of complex diseases but also enable the screening of individuals at high risk of complex diseases,  providing opportunities for prophylactic actions, and guiding precision medicine.  In healthcare practice, however, there are still huge gaps to be bridged between genetic discoveries and their actual uses. One major obstacle to filling the gap lies in the translation of individual genetic profiles into actionable healthcare decisions. In most cases, associations between genetic variants and phenotypes are not as simple as those in exceptionally rare monogenic disorders, which are genetically controlled by a single locus (e.g., a single base‐pair substitution).  Most disorders are determined by the overall effect of numerous genetic variations with varying degrees of contributions.  This complex genetic architecture makes it clinically challenging to translate implicit genetic information into explicit healthcare decisions.  Current genotyping techniques acquire genetic data by specialized instrumentation and data analysis pipeline,  but of particular note, it is the phenotype predicted from multiplex genetic data that offers health benefits. In the context of precision medicine, for example, clinical prescribing actions rely on individual drug responses predicted from multiple genetic determinants of drug efficacy and/or toxicity.  Accordingly, in silico computational tools,  along with documented guidelines  and expert data interpretation,  are quite necessary to finalize the workflow of traditional genetic services but are inconsistently implemented across clinical laboratories,  thereby discouraging wide‐scale incorporation of genetics‐directed personalized healthcare into routine clinical practice. A theoretically promising approach for intelligent genotype‐phenotype interpretation is DNA‐based molecular computing circuitry. Thus far, DNA has been extensively used in molecular engineering as a versatile, programmable, and biocompatible building block.  Furthermore, the thermodynamics and kinetics of DNA hybridization and strand displacement have been thoroughly characterized,  making it possible to tailor DNA circuits for diverse computational tasks, such as Boolean logic operations,  basic arithmetic computation,  and even neural network mimicry.  More importantly, DNA is much flexible to interface with biological systems, opening up many possibilities for DNA computation in clinical diagnostics,  bio‐imaging,  and therapeutics.  Therefore, it can be anticipated to engineer an intelligent DNA decoding circuit able to recognize and interpret multiplex genetic profiles into actionable healthcare decisions in a “tube” without in silico computation and expert data interpretation. Herein, we established a general framework built on a programmable DNA decoding circuit to facilitate the translation of complex genetic profiles into actionable medical decisions. To make genetic information recognizable, we developed a ligase‐based signal transforming strategy to accurately and conveniently convert multiplex allelic information into a battery of single‐stranded DNA uniquely mapped to DNA gates in the DNA decoding circuit by a one‐pot ligase‐dependent reaction (LDR). This strategy leverages the allele‐discriminating power of thermostable ligase, thus circumventing mis‐extension errors often encountered in polymerase‐based genotyping methods, such as TaqMan genotyping assays and amplification refractory mutation system polymerase chain reaction (PCR).  The DNA decoding circuit is a highly modular and programmable arithmetic circuit with the ability to “acquire” genetic knowledge and then translate various sets of genetic binary codes (i.e., the presence of an allele is coded as 1 and absence as 0) into actionable healthcare decisions (indicated by fluorescent signals) according to pre‐defined programs with modest instrument requirements. We show that the DNA decoding circuit implements expected computation for varying combinations of synthetic inputs and when tested in human genomic DNA samples, still achieves desirable performance on the interpretation of “genetic codes” of CYP2C19 , a pharmacogenetics (PGx) gene, into drug responses, with accuracy equivalent to that of Sanger sequencing. We also demonstrate that circuit layouts can be readily rewired for other healthcare purposes, such as interpreting pharmacogenetics tests, providing drug dosing recommendations, and implementing reliable molecular calculation of polygenic risk score (PRS) as well as PRS‐informed disease risk assessment.   Results 2.1 Overview of the DNA‐Powered Intelligent Genetic Decoder A DNA‐based decoding circuit with previously acquired genetics knowledge can recognize multiplex genetic profiles and then output the most probable phenotype of an individual in a tube, thus it can provide users with explicit reports in a “sample in, answer out” fashion ( Figure ). Also, circuit programs are reconfigurable for varying healthcare purposes, such as dose recommendations and disease risk assessment. The underlying principle is a series of molecular arithmetic operations over a battery of effect alleles with established roles in inter‐individual phenotypic variabilities (Figure ). At first, all alleles are multiplied individually by their respective weights (Figure ), which depend on the contribution to phenotype expression. Specifically, a loss‐of‐function allele may take a negative weight whilst a gain‐of‐function allele may take a positive weight. All weighted alleles then add up to an overall score (Figure ), which varies with distinct phenotypic classes. Accordingly, comparators are incorporated to distinguish different levels of weighted‐sum and output fluorescent signals (Figure ). To program a decoding circuit, a critical step is to learn genetic features underlying diverse phenotypes, that is, to look for an optimal array of multiplicative coefficients assigned to each allele. An optimal coefficient array should make all output classes separable in the weighted‐sum space such that the classification task is formulated as a set of mathematical operations. To clarify it, we take four phenotypic classes of drug responses caused by CYP2C19 single nucleotide polymorphisms (SNPs) as a demonstration.  Since each SNP (or exactly biallelic SNP) gives rise to two alternative alleles within a specific locus, multiple SNPs across different loci can lead to various combinations of all alternative alleles, which we refer to as “genetic codes” (Figure ). All genotypes were denoted by the star  alleles according to the star allele nomenclature, where each number represents a specific allele with altered functionality caused by a single variant or a group of genetic variations.  For each CYP2C19 genotype, there is a particular set of “genetic codes” encoded by six possible CYP2C19 alleles (Table , Supporting Information). For functional classification, individuals with different genotypes are generally categorized into four phenotypic categories: poor (PM), intermediate (IM), extensive (EM), and ultrarapid (UM) metabolizers. Genotypes of the same phenotypic class, though holding distinct “genetic codes”, should well separate from other classes in the 2D weighted‐sum space (Figure ). Accordingly, multiplex genetic profiles were reduced to two dimensions ( X vs Y ), and the classification of multiple phenotypic classes thus can be achieved by comparing weighted‐sum numerically with reference values (Figure ). To facilitate the learning process of allelic features, a MATLAB script was written to automatically look for fitted coefficient arrays for a given phenotype‐genotype dataset (Note , Supporting Information). 2.2 Conversion from Genetic Profiles into Recognizable Inputs Converting genetic profiles into recognizable inputs is an essential prerequisite in attempts to engineer molecular genetic decoders. To this end, we developed a thermostable ligase‐based transforming strategy to accurately transform multiplex genetic information into DNA strands with unique sequences through a one‐pot ligase‐dependent reaction (LDR) ( Figure ). LDR achieves precise allelic discrimination by the thermostable 9°N™ ligase that can covalently join the end of two immediately adjacent and perfectly matched primers following binding to the target template at elevated temperatures, whereas any variations (e.g., SNPs, Insertions, and Deletions) at the junction can prevent them from being successfully ligated (Figure ). Each allele has an upstream allele‐specific primer with a particular toehold domain at its 5’ end, and a phosphorylated downstream primer (Figure ). When the target allele is present, both primers will hybridize with the target allele and be joined together by ligase, generating a linked primer at each thermal cycle (Figure ). However, there is none of linked primers in the absence of the target allele (Figure ). The allele‐discriminating power of LDR was cross‐validated using five primer pairs designed for five CYP2C19 alleles (see Table , Supporting Information for detailed sequences). Five synthetic plasmids carrying a 4.81 kb gene fragment of corresponding alleles were used as the template of LDR. It demonstrated that linked primers were detectable only in the context of targeted allelic templates, whereas none of linked primers were detected in the context of non‐targeted allelic templates (Figure ). Next, we used ImageJ software to quantify band intensities of linked and unlinked primers on the gel image and approximated their concentrations based on the initial total concentration,  revealing that the ligation efficiency ranged from 35% to 64% and unspecific ligation events were no more than 2.1% (Table , Supporting Information). Overall, it demonstrates that LDR maintains a high selectivity against point mutations while giving a desirable yield of ligated products. As the toehold domain on the allele‐specific primer is perfectly complementary to that of an allele‐specific DNA arithmetic gate, successfully linked allele‐specific primers can initiate a strand displacement reaction with allele‐specific DNA gate and displace an output strand (Figure ). Instead, both unlinked primers are thermodynamically unfavorable to displacing the output (Figure ). Experimentally characterized kinetics behaviors of five LDR primer pairs verified that only linked primer can initiate the strand displacement reaction, yet unlinked primers were observed to be unable to complete the strand displacement (Figure ) (Figure , Supporting Information). Characterization of reaction products by gel electrophoresis also shows that only linked primer is able to displace the output strand from the DNA gate and hybridize with the bottom strand of the DNA gate (Figure , Supporting Information), further verifying this proposed reaction pathway. We also demonstrated that the ligase‐based method could reliably transform DNA methylation into unique circular ssDNA products when combined with molecular inversion probe (MIP). MIP is essentially a single‐stranded oligonucleotide, with two ends adjacent and complementary to flanking sequences of target CpG, and a middle domain for downstream mathematical operations. Both ends of the probe targeting the methylation sites is complementary to the methylated template converted by sodium sulfite (Figure , Supporting Information). We designed a pair of MIPs for Cg16655791 and experimentally characterized that linear MIPs will be circularized by thermostable ligase only in the context of correct template, while mismatches at the joint will prevent it from being ligated. Exonuclease I and exonuclease III will remove linear probes and double‐stranded template respectively, leaving only circular probes for downstream molecular operations (Figure , Supporting Information). 2.3 Molecular Arithmetic Operations over Genetic Profiles The general scheme for molecular mathematic implementation is based on the toehold‐mediated strand‐displacement reaction network between distinct DNA arithmetic gates at predetermined concentrations. Specifically, multiplicative computation was powered by a panel of “seesaw catalyst”‐based multiplication gates designated for each allele (Figure , Supporting Information). Though concentrations of linked primers may vary with different alleles, weighted outputs for each allele can be definitive by fixing the concentration of the multiplication gate in that a small number of input species can activate catalytic strand‐displacement reaction.  The multiplication gate has a double‐stranded central domain and two single‐stranded flanking domains, of which one side is the toehold that can be recognized by a linked primer to initiate branch migration, while the other shares identical sequences with the branch migration domain of downstream summation gates ( Figure ). After hybridizing with the toehold on the multiplication gate, the linked primer will displace an output strand and activate the other toehold (d*) simultaneously (Figure ). Excessive fuel strands will recognize the newly uncovered toehold and initiate another strand‐displacement cycle to free up the linked primer, thereby allowing linked primers to be reused over multiple cycles (Figure ). Outputs of DNA multiplication gates at user‐defined molar concentrations were demonstrated to be completely displaced by corresponding linked primers even at six‐times lower quantity (Figure ) (Figure , Supporting Information). The summation of weighted species displaced from multiplication gates was realized by two partially double‐stranded DNA summation gates, which respectively convert all positively weighted species into identical output X and negatively‐weighted species into Y (Figure ). With a specific toehold domain, any weighted species with complementary sequences can initiate branch migration on summation gate (Figure ). Each summation gate was individually tested with synthetic inputs of multiplication gates, demonstrating that the summation gate could exactly transform weighted species displaced from distinct upstream multiplication gates into identical outputs (Figure ). We next validated whether the two arithmetic gates could accurately compute weighted‐sum for different CYP2C19 genotypes. With a special panel of “genetic codes,” each genotype would have a particular combination of linked primers in LDR products. Using serial combinations of synthetic linked primers as inputs, the circuit yielded desired concentrations of weighted‐sum species X and Y for all test cases (Figure , Supporting Information), indicating that allele‐specific multiplication gates could work in parallel to compute weighted species. When input species were changed from synthetic DNA strands to products of LDR (plasmid templates were prepared according to the specific allele status of each genotype, as shown in Table , Supporting Information), the circuit still outputted desired concentrations of weighted‐sum X and Y for all tested CYP2C19 genotypes (Figure ) (Figure , Supporting Information), indicating that the circuit could recognize transformed “genetic codes” and perform arithmetic operations over the mixture of ligase‐mediated linked primers, whose components may vary with genotypes. However, it still requires a subtraction gate to compute “ X+Y ” in that species Y represents a negative weighted‐sum. This can be achieved by a cooperative hybridization scheme,  whereby stoichiometrically equivalent amounts of species X and Y can simultaneously hybridize with the subtraction gate and irreversibly annihilate each other only if both species are present. To realize molecular subtractive operations, we additionally designed two DNA transformation gates (Figure , Supporting Information), whose strand‐displacement kinetics were programmed to be much slower than the subtraction gate by cutting 2 nucleotides from its toehold domains (the length and sequence of the toehold domain are known to significantly influence the kinetics of strand displacement  ). Accordingly, the weighted‐sum with a lower absolute quantity will be exhausted and the leftover of the higher weighted‐sum species will sequentially hybridize with and displace an output strand (H1 or H2) from a transformation gate stoichiometrically (Figure , Supporting Information). We used 12 combinations of weighted‐sum X and Y at initial concentrations ranging from 0 to 120 n m to demonstrate the subtraction scheme and observed that the circuit behavior generally agree with the proposed design: if the initial concentration of Y is higher than that of X , then the concentration of the H1 strand is approximate to the leftover amount of Y ([H1] ∞ ≈[ Y ] 0 −[ X ] 0 ); when the initial concentration of X is higher that of Y , the concentration of H2 strand is approximate to the leftover amount of X ([H2] ∞ ≈[ X ] 0 −[ Y ] 0 ) (Figure , Supporting Information). In this way, weighted‐sum species computed by upstream arithmetic gates are transformed into different concentrations of species H1 or H2 that can be processed by downstream comparators. To implement numeric comparison with DNA molecules, we took advantage of two threshold gates in conjunction with four catalytic fluorescence reporting systems, among which VIC and Cy5 reporters were respectively programmed to indicate the presence of H1 (H1 > 0) and H2 (H2 > 0), while FAM and ROX reporters were expected to output fluorescence only if H1 and H2 exceed a user‐defined threshold, respectively (Figure , Supporting Information). The species H1 and H2 were programmed to preferentially hybridize with threshold gates by supplementing the toehold on threshold gates with extra nucleotides, thereby rendering strand displacement reactions of threshold gates much faster than those of restoration gates (Figure , Supporting Information). We experimentally tested the thresholding effectiveness using varying concentrations of species H and observed that the circuit could trigger the desired fluorescent signal(s) for species H above the threshold, whereas undesired signals, though detectable, maintained a desirable OFF state in cases of species H below the threshold (Figure , Supporting Information). Further, we verified the thresholding effectiveness of the comparator connected with a subtraction module using six distinct combinations of X and Y (Figure ). In agreement with the circuit design, significant VIC signals were only observed in instances of [ X ] 0 −[ Y ] 0 < 0 ([H1] > 0), whilst Cy5 signals were observed only if [ X ] 0 −[ Y ] 0 > 0 ([H2] > 0), whereas ROX and FAM signals remarkably increased only when [ X ] 0 −[ Y ] 0 exceeds the predefined threshold (Figure ). Overall, these results suggest that the DNA‐based decoding circuit can recognize and compute genetic profiles of distinct genotypes and export desired fluorescent signals predefined for the phenotypic class of a genotype. 2.4 Evaluating the Performance of the Genetic Decoder on Human Genomic DNA Samples To demonstrate the robustness and usefulness of the genetic decoder, we collected 30 distinct human genomic DNA samples whose CYP2C19 genotypes have been validated using Sanger sequencing (Table , Supporting Information). Initially, we used synthetic inputs to test if the decoding circuit could yield expected fluorescent patterns for each phenotype. Though desired fluorescent patterns were observed for all phenotypic classes, undesired signals, particularly those of FAM and Cy5, remained detectable (Figure , Supporting Information). To reproducibly demarcate ON/OFF state, an RFU threshold of 0.4 was defined according to FAM/Cy5 RFU of simulated intermediate metabolizer samples. Further, we evaluated the performance of the entire workflow using binary combinations of synthetic plasmids as simulated samples. As observed, all genotypes could be assigned to the expected phenotype based on their fluorescent outputs ( Figure ). Finally, we assessed the performance of the DNA‐powered genetic decoder on 30 human genomic DNA samples. All samples were first subjected to pre‐amplification owing to the extremely low copy number. Premixed circuit components were then pooled into the finished LDR. Without any post‐reaction treatment, all testing procedures could be accomplished within a test tube automatically. As shown in the experimental data, genotypes of the same phenotypic class triggered a similar fluorescence pattern within 2.5 h, which agrees with the expected fluorescence pattern programmed for each class, more importantly, samples tested in different batches could also yield expected fluorescence patterns (Figure ). These results suggest the robustness and usefulness of the DNA genetic decoder in real genomic DNA samples. Moreover, we have experimentally demonstrated that the ligase‐based method could transform rare nucleic acid variants with low allele frequencies effectively when combined with blocker displacement amplification (Figure , Supporting Information).  Rare alleles with even 1% VAF were demonstrated to be effectively enriched through BDA (Figure , Supporting Information), and can trigger accurate multiplicative operations after signal transformation via LDR (Figure , Supporting Information) 2.5 Reprogramming DNA Decoding Circuits for Diverse Healthcare Purposes After experimentally verifying the effectiveness of the DNA‐powered decoding circuit engineered for decoding CYP2C19 “genetic codes” of human genomic DNA samples, we wondered whether the DNA‐based genetic decoder could be readily reprogrammed for diverse healthcare purposes. Since all components are simple DNA duplexes and highly modular, we expect that programs for more complicated tasks can be reconfigured by simply rewiring or tuning concentrations of circuit modules. In our initial design, the decoding circuit was programmed to interpret 4 phenotypic classes from various sets of CYP2C19 “genetic codes.” For other PGx‐genes, phenotypic classes may be annotated in different ways. In the case of DPYD , there are three likely phenotypic classes (normal, intermediate, poor), which are classified based on activity scores of DPYD genotypes.  Using the DPYD genotype‐phenotype dataset as input, an optimal array of allelic coefficients can be identified by in silico learning. With the fitted coefficient array, three phenotypic classes were observed to be well separated in the weighted‐sum space ( Figure ). From this, we reprogrammed a decoding circuit to interpret 30 sets of DPYD “genetic codes” into 3 phenotypic classes based on 14 alleles caused by 7 DPYD SNPs (Table , Supporting Information) and then predicted circuit behaviors using the CRN simulator.  The simulation can model multiple sets of strand‐displacement reactions and predict circuit behaviors to facilitate circuit optimization.  The circuit was predicted to perform desired computation for the allele‐specific input set of each genotype and generate the correct fluorescence pattern programmed for the phenotypic class of each genotype (Figure , Supporting Information). When synthetic input combinations of DPYD genotypes were introduced, the circuit was experimentally observed to output expected fluorescence patterns (Figure ). Furthermore, we expected that the circuit output could shift from the phenotype of drug responses to the dose of a specific medication. Phenotypic variations caused by a single PGx gene may affect the usage of numerous drugs, and likewise, the dose of a single drug is often influenced by several PGx genes. Therefore, we explored the potential to provide dosage recommendations on thiopurines, whose effectiveness and toxicity have been verified to interact with two PGx genes— TPMT and NUDT15 .  To reprogram a decoding circuit for dose recommendations, we identified coefficient arrays for 8 TPMT and 4 NUDT15 alleles individually using the abovementioned script (Table , Supporting Information). The circuit was designed to calculate a weighted‐sum for genotypes of both genes in parallel, compare the weighted‐sum separately when either weighted‐sum is negative (X1+Y1<0 or X2+Y2<0) but compare the weighted‐sum collectively if both weighted‐sum are positive (X1+Y1>0 and X2+Y2>0) (Figure ). This is because the normal dose is only applicable to normal metabolizers (X1+Y1>0 and X2+Y2>0) of TPMT and NUDT15, whereas dose reductions are recommended for poor and intermediate metabolizers (X1+Y1>0 or X2+Y2>0) of either TPMT or NUDT15 (Table , Supporting Information).  Circuit simulations predicted correct fluorescent patterns for input combinations of all given genotypes (Figure , Supporting Information), and experimental validation using synthetic input combinations of given DPYD genotypes demonstrates expected fluorescence patterns programmed for each recommended dosage (Figure ). Another two linear models were also successfully trained for interpreting SLCO1B1 and UGT1A1 genetic profiles into dose recommendations (Figure , Supporting Information), suggesting that the DNA decoding circuit is highly programmable to guide personalized healthcare. In addition to guiding drug regimes, another rapidly emerging precision medicine tool is the polygenic risk score (PRS), which estimates a person's genetic risk of developing a particular disease by cumulating the effects of various risk alleles with relatively small individual effects on incidence.  Though the clinical validity and utility of polygenic risk profiling for risk stratification, screening, and prediction have been evaluated in detail,  it remains challenging to implement PRS in routine clinical and public health practice,  which prompted us to explore the utility of DNA decoding circuit in the molecular calculation of PRS and PRS‐informed disease screening, considering its high programmability and modularized circuit layout. We thus abstracted a PGS scoring system for pancreatic cancer from the Polygenic Score Catalog ( http://www.pgscatalog.org ).  This PGS scoring system includes 22 risk variants with varying effect sizes (Odds Ratio, OR), and PRS could be calculated as a weighted‐sum score of all effect alleles using their effect sizes from GWAS as weights ( Figure ). Notably, the allelic dosage was defined by the allele copy number (that is, “1” for a heterozygote and “2” for a homozygote), whereas each allele has a binary state of 0 or 1 in our previously engineered circuit. To achieve a ternary input mode, the multiplicative coefficient of effect alleles was programmed to double its original lnOR while the coefficient of other alleles was set to have a negative lnOR (Figure ). In this way, the concentration of weighted species is 2lnOR if the allele is homozygous but reduces to lnOR if the allele is heterozygous, and there are none of positively weighted species in the absence of the allele. Weighted species of each risk locus were then integrated by a summation gate, thereby setting the stage for the molecule calculation of PRS (Figure ). Using synthetic input patterns to simulate homozygous/heterozygous allele samples, we experimentally demonstrated that the DNA circuit could implement reliable molecular calculation of risk score for individual risk locus (Figure , Supporting Information). We also simulated the circuit behavior with the genotypic dataset of 2504 individuals from the 1000 Genomes Project Phase 3 and statistically compared output concentrations of weighted‐sum species with their supposed PRS(s), which were observed to be normally distributed (mean: 4.30; standard deviation: 0.69) (Figure ). Risk stratification based on supposed PRS(s) of the 2504 individuals resulted in three categories: low (bottom decile), intermediate (2nd–8th decile), and high (top decile) (Figure ), and we noted that the stoichiometric quantity of the weighted‐sum output was near to the supposed PRS for most test cases and demonstrated a similar distribution with that of supposed PRS (mean: 4.27; standard deviation: 0.34) (Figure ). Accordingly, we reasoned that high‐risk individual with top10% PRS could be efficiently screened by fluorescent outputs when a comparator module was incorporated into the circuit (Figure ). Using synthetic input patterns of individuals with varying PRS levels (Table , Supporting Information), circuit behaviors were expectedly demonstrated to have an ON trajectory in test cases with top1% PRS and top5% PRS, while maintaining an OFF state for test cases with intermediate PRS (Figure ), suggesting that the DNA‐powered genetic decoder has great potential to efficiently and straightforwardly implement PRS‐informed individualized disease risk prediction and perform populational risk stratification. Compared with the well‐established Sanger sequencing for genetic testing and interpretation (Table , Supporting Information), our DNA‐computation‐based platform shows advantageous features including high multiplexity, automatability, and target expandability, as well as fast turnaround time.   Overview of the DNA‐Powered Intelligent Genetic Decoder A DNA‐based decoding circuit with previously acquired genetics knowledge can recognize multiplex genetic profiles and then output the most probable phenotype of an individual in a tube, thus it can provide users with explicit reports in a “sample in, answer out” fashion ( Figure ). Also, circuit programs are reconfigurable for varying healthcare purposes, such as dose recommendations and disease risk assessment. The underlying principle is a series of molecular arithmetic operations over a battery of effect alleles with established roles in inter‐individual phenotypic variabilities (Figure ). At first, all alleles are multiplied individually by their respective weights (Figure ), which depend on the contribution to phenotype expression. Specifically, a loss‐of‐function allele may take a negative weight whilst a gain‐of‐function allele may take a positive weight. All weighted alleles then add up to an overall score (Figure ), which varies with distinct phenotypic classes. Accordingly, comparators are incorporated to distinguish different levels of weighted‐sum and output fluorescent signals (Figure ). To program a decoding circuit, a critical step is to learn genetic features underlying diverse phenotypes, that is, to look for an optimal array of multiplicative coefficients assigned to each allele. An optimal coefficient array should make all output classes separable in the weighted‐sum space such that the classification task is formulated as a set of mathematical operations. To clarify it, we take four phenotypic classes of drug responses caused by CYP2C19 single nucleotide polymorphisms (SNPs) as a demonstration.  Since each SNP (or exactly biallelic SNP) gives rise to two alternative alleles within a specific locus, multiple SNPs across different loci can lead to various combinations of all alternative alleles, which we refer to as “genetic codes” (Figure ). All genotypes were denoted by the star  alleles according to the star allele nomenclature, where each number represents a specific allele with altered functionality caused by a single variant or a group of genetic variations.  For each CYP2C19 genotype, there is a particular set of “genetic codes” encoded by six possible CYP2C19 alleles (Table , Supporting Information). For functional classification, individuals with different genotypes are generally categorized into four phenotypic categories: poor (PM), intermediate (IM), extensive (EM), and ultrarapid (UM) metabolizers. Genotypes of the same phenotypic class, though holding distinct “genetic codes”, should well separate from other classes in the 2D weighted‐sum space (Figure ). Accordingly, multiplex genetic profiles were reduced to two dimensions ( X vs Y ), and the classification of multiple phenotypic classes thus can be achieved by comparing weighted‐sum numerically with reference values (Figure ). To facilitate the learning process of allelic features, a MATLAB script was written to automatically look for fitted coefficient arrays for a given phenotype‐genotype dataset (Note , Supporting Information).   Conversion from Genetic Profiles into Recognizable Inputs Converting genetic profiles into recognizable inputs is an essential prerequisite in attempts to engineer molecular genetic decoders. To this end, we developed a thermostable ligase‐based transforming strategy to accurately transform multiplex genetic information into DNA strands with unique sequences through a one‐pot ligase‐dependent reaction (LDR) ( Figure ). LDR achieves precise allelic discrimination by the thermostable 9°N™ ligase that can covalently join the end of two immediately adjacent and perfectly matched primers following binding to the target template at elevated temperatures, whereas any variations (e.g., SNPs, Insertions, and Deletions) at the junction can prevent them from being successfully ligated (Figure ). Each allele has an upstream allele‐specific primer with a particular toehold domain at its 5’ end, and a phosphorylated downstream primer (Figure ). When the target allele is present, both primers will hybridize with the target allele and be joined together by ligase, generating a linked primer at each thermal cycle (Figure ). However, there is none of linked primers in the absence of the target allele (Figure ). The allele‐discriminating power of LDR was cross‐validated using five primer pairs designed for five CYP2C19 alleles (see Table , Supporting Information for detailed sequences). Five synthetic plasmids carrying a 4.81 kb gene fragment of corresponding alleles were used as the template of LDR. It demonstrated that linked primers were detectable only in the context of targeted allelic templates, whereas none of linked primers were detected in the context of non‐targeted allelic templates (Figure ). Next, we used ImageJ software to quantify band intensities of linked and unlinked primers on the gel image and approximated their concentrations based on the initial total concentration,  revealing that the ligation efficiency ranged from 35% to 64% and unspecific ligation events were no more than 2.1% (Table , Supporting Information). Overall, it demonstrates that LDR maintains a high selectivity against point mutations while giving a desirable yield of ligated products. As the toehold domain on the allele‐specific primer is perfectly complementary to that of an allele‐specific DNA arithmetic gate, successfully linked allele‐specific primers can initiate a strand displacement reaction with allele‐specific DNA gate and displace an output strand (Figure ). Instead, both unlinked primers are thermodynamically unfavorable to displacing the output (Figure ). Experimentally characterized kinetics behaviors of five LDR primer pairs verified that only linked primer can initiate the strand displacement reaction, yet unlinked primers were observed to be unable to complete the strand displacement (Figure ) (Figure , Supporting Information). Characterization of reaction products by gel electrophoresis also shows that only linked primer is able to displace the output strand from the DNA gate and hybridize with the bottom strand of the DNA gate (Figure , Supporting Information), further verifying this proposed reaction pathway. We also demonstrated that the ligase‐based method could reliably transform DNA methylation into unique circular ssDNA products when combined with molecular inversion probe (MIP). MIP is essentially a single‐stranded oligonucleotide, with two ends adjacent and complementary to flanking sequences of target CpG, and a middle domain for downstream mathematical operations. Both ends of the probe targeting the methylation sites is complementary to the methylated template converted by sodium sulfite (Figure , Supporting Information). We designed a pair of MIPs for Cg16655791 and experimentally characterized that linear MIPs will be circularized by thermostable ligase only in the context of correct template, while mismatches at the joint will prevent it from being ligated. Exonuclease I and exonuclease III will remove linear probes and double‐stranded template respectively, leaving only circular probes for downstream molecular operations (Figure , Supporting Information).   Molecular Arithmetic Operations over Genetic Profiles The general scheme for molecular mathematic implementation is based on the toehold‐mediated strand‐displacement reaction network between distinct DNA arithmetic gates at predetermined concentrations. Specifically, multiplicative computation was powered by a panel of “seesaw catalyst”‐based multiplication gates designated for each allele (Figure , Supporting Information). Though concentrations of linked primers may vary with different alleles, weighted outputs for each allele can be definitive by fixing the concentration of the multiplication gate in that a small number of input species can activate catalytic strand‐displacement reaction.  The multiplication gate has a double‐stranded central domain and two single‐stranded flanking domains, of which one side is the toehold that can be recognized by a linked primer to initiate branch migration, while the other shares identical sequences with the branch migration domain of downstream summation gates ( Figure ). After hybridizing with the toehold on the multiplication gate, the linked primer will displace an output strand and activate the other toehold (d*) simultaneously (Figure ). Excessive fuel strands will recognize the newly uncovered toehold and initiate another strand‐displacement cycle to free up the linked primer, thereby allowing linked primers to be reused over multiple cycles (Figure ). Outputs of DNA multiplication gates at user‐defined molar concentrations were demonstrated to be completely displaced by corresponding linked primers even at six‐times lower quantity (Figure ) (Figure , Supporting Information). The summation of weighted species displaced from multiplication gates was realized by two partially double‐stranded DNA summation gates, which respectively convert all positively weighted species into identical output X and negatively‐weighted species into Y (Figure ). With a specific toehold domain, any weighted species with complementary sequences can initiate branch migration on summation gate (Figure ). Each summation gate was individually tested with synthetic inputs of multiplication gates, demonstrating that the summation gate could exactly transform weighted species displaced from distinct upstream multiplication gates into identical outputs (Figure ). We next validated whether the two arithmetic gates could accurately compute weighted‐sum for different CYP2C19 genotypes. With a special panel of “genetic codes,” each genotype would have a particular combination of linked primers in LDR products. Using serial combinations of synthetic linked primers as inputs, the circuit yielded desired concentrations of weighted‐sum species X and Y for all test cases (Figure , Supporting Information), indicating that allele‐specific multiplication gates could work in parallel to compute weighted species. When input species were changed from synthetic DNA strands to products of LDR (plasmid templates were prepared according to the specific allele status of each genotype, as shown in Table , Supporting Information), the circuit still outputted desired concentrations of weighted‐sum X and Y for all tested CYP2C19 genotypes (Figure ) (Figure , Supporting Information), indicating that the circuit could recognize transformed “genetic codes” and perform arithmetic operations over the mixture of ligase‐mediated linked primers, whose components may vary with genotypes. However, it still requires a subtraction gate to compute “ X+Y ” in that species Y represents a negative weighted‐sum. This can be achieved by a cooperative hybridization scheme,  whereby stoichiometrically equivalent amounts of species X and Y can simultaneously hybridize with the subtraction gate and irreversibly annihilate each other only if both species are present. To realize molecular subtractive operations, we additionally designed two DNA transformation gates (Figure , Supporting Information), whose strand‐displacement kinetics were programmed to be much slower than the subtraction gate by cutting 2 nucleotides from its toehold domains (the length and sequence of the toehold domain are known to significantly influence the kinetics of strand displacement  ). Accordingly, the weighted‐sum with a lower absolute quantity will be exhausted and the leftover of the higher weighted‐sum species will sequentially hybridize with and displace an output strand (H1 or H2) from a transformation gate stoichiometrically (Figure , Supporting Information). We used 12 combinations of weighted‐sum X and Y at initial concentrations ranging from 0 to 120 n m to demonstrate the subtraction scheme and observed that the circuit behavior generally agree with the proposed design: if the initial concentration of Y is higher than that of X , then the concentration of the H1 strand is approximate to the leftover amount of Y ([H1] ∞ ≈[ Y ] 0 −[ X ] 0 ); when the initial concentration of X is higher that of Y , the concentration of H2 strand is approximate to the leftover amount of X ([H2] ∞ ≈[ X ] 0 −[ Y ] 0 ) (Figure , Supporting Information). In this way, weighted‐sum species computed by upstream arithmetic gates are transformed into different concentrations of species H1 or H2 that can be processed by downstream comparators. To implement numeric comparison with DNA molecules, we took advantage of two threshold gates in conjunction with four catalytic fluorescence reporting systems, among which VIC and Cy5 reporters were respectively programmed to indicate the presence of H1 (H1 > 0) and H2 (H2 > 0), while FAM and ROX reporters were expected to output fluorescence only if H1 and H2 exceed a user‐defined threshold, respectively (Figure , Supporting Information). The species H1 and H2 were programmed to preferentially hybridize with threshold gates by supplementing the toehold on threshold gates with extra nucleotides, thereby rendering strand displacement reactions of threshold gates much faster than those of restoration gates (Figure , Supporting Information). We experimentally tested the thresholding effectiveness using varying concentrations of species H and observed that the circuit could trigger the desired fluorescent signal(s) for species H above the threshold, whereas undesired signals, though detectable, maintained a desirable OFF state in cases of species H below the threshold (Figure , Supporting Information). Further, we verified the thresholding effectiveness of the comparator connected with a subtraction module using six distinct combinations of X and Y (Figure ). In agreement with the circuit design, significant VIC signals were only observed in instances of [ X ] 0 −[ Y ] 0 < 0 ([H1] > 0), whilst Cy5 signals were observed only if [ X ] 0 −[ Y ] 0 > 0 ([H2] > 0), whereas ROX and FAM signals remarkably increased only when [ X ] 0 −[ Y ] 0 exceeds the predefined threshold (Figure ). Overall, these results suggest that the DNA‐based decoding circuit can recognize and compute genetic profiles of distinct genotypes and export desired fluorescent signals predefined for the phenotypic class of a genotype.   Evaluating the Performance of the Genetic Decoder on Human Genomic DNA Samples To demonstrate the robustness and usefulness of the genetic decoder, we collected 30 distinct human genomic DNA samples whose CYP2C19 genotypes have been validated using Sanger sequencing (Table , Supporting Information). Initially, we used synthetic inputs to test if the decoding circuit could yield expected fluorescent patterns for each phenotype. Though desired fluorescent patterns were observed for all phenotypic classes, undesired signals, particularly those of FAM and Cy5, remained detectable (Figure , Supporting Information). To reproducibly demarcate ON/OFF state, an RFU threshold of 0.4 was defined according to FAM/Cy5 RFU of simulated intermediate metabolizer samples. Further, we evaluated the performance of the entire workflow using binary combinations of synthetic plasmids as simulated samples. As observed, all genotypes could be assigned to the expected phenotype based on their fluorescent outputs ( Figure ). Finally, we assessed the performance of the DNA‐powered genetic decoder on 30 human genomic DNA samples. All samples were first subjected to pre‐amplification owing to the extremely low copy number. Premixed circuit components were then pooled into the finished LDR. Without any post‐reaction treatment, all testing procedures could be accomplished within a test tube automatically. As shown in the experimental data, genotypes of the same phenotypic class triggered a similar fluorescence pattern within 2.5 h, which agrees with the expected fluorescence pattern programmed for each class, more importantly, samples tested in different batches could also yield expected fluorescence patterns (Figure ). These results suggest the robustness and usefulness of the DNA genetic decoder in real genomic DNA samples. Moreover, we have experimentally demonstrated that the ligase‐based method could transform rare nucleic acid variants with low allele frequencies effectively when combined with blocker displacement amplification (Figure , Supporting Information).  Rare alleles with even 1% VAF were demonstrated to be effectively enriched through BDA (Figure , Supporting Information), and can trigger accurate multiplicative operations after signal transformation via LDR (Figure , Supporting Information)   Reprogramming DNA Decoding Circuits for Diverse Healthcare Purposes After experimentally verifying the effectiveness of the DNA‐powered decoding circuit engineered for decoding CYP2C19 “genetic codes” of human genomic DNA samples, we wondered whether the DNA‐based genetic decoder could be readily reprogrammed for diverse healthcare purposes. Since all components are simple DNA duplexes and highly modular, we expect that programs for more complicated tasks can be reconfigured by simply rewiring or tuning concentrations of circuit modules. In our initial design, the decoding circuit was programmed to interpret 4 phenotypic classes from various sets of CYP2C19 “genetic codes.” For other PGx‐genes, phenotypic classes may be annotated in different ways. In the case of DPYD , there are three likely phenotypic classes (normal, intermediate, poor), which are classified based on activity scores of DPYD genotypes.  Using the DPYD genotype‐phenotype dataset as input, an optimal array of allelic coefficients can be identified by in silico learning. With the fitted coefficient array, three phenotypic classes were observed to be well separated in the weighted‐sum space ( Figure ). From this, we reprogrammed a decoding circuit to interpret 30 sets of DPYD “genetic codes” into 3 phenotypic classes based on 14 alleles caused by 7 DPYD SNPs (Table , Supporting Information) and then predicted circuit behaviors using the CRN simulator.  The simulation can model multiple sets of strand‐displacement reactions and predict circuit behaviors to facilitate circuit optimization.  The circuit was predicted to perform desired computation for the allele‐specific input set of each genotype and generate the correct fluorescence pattern programmed for the phenotypic class of each genotype (Figure , Supporting Information). When synthetic input combinations of DPYD genotypes were introduced, the circuit was experimentally observed to output expected fluorescence patterns (Figure ). Furthermore, we expected that the circuit output could shift from the phenotype of drug responses to the dose of a specific medication. Phenotypic variations caused by a single PGx gene may affect the usage of numerous drugs, and likewise, the dose of a single drug is often influenced by several PGx genes. Therefore, we explored the potential to provide dosage recommendations on thiopurines, whose effectiveness and toxicity have been verified to interact with two PGx genes— TPMT and NUDT15 .  To reprogram a decoding circuit for dose recommendations, we identified coefficient arrays for 8 TPMT and 4 NUDT15 alleles individually using the abovementioned script (Table , Supporting Information). The circuit was designed to calculate a weighted‐sum for genotypes of both genes in parallel, compare the weighted‐sum separately when either weighted‐sum is negative (X1+Y1<0 or X2+Y2<0) but compare the weighted‐sum collectively if both weighted‐sum are positive (X1+Y1>0 and X2+Y2>0) (Figure ). This is because the normal dose is only applicable to normal metabolizers (X1+Y1>0 and X2+Y2>0) of TPMT and NUDT15, whereas dose reductions are recommended for poor and intermediate metabolizers (X1+Y1>0 or X2+Y2>0) of either TPMT or NUDT15 (Table , Supporting Information).  Circuit simulations predicted correct fluorescent patterns for input combinations of all given genotypes (Figure , Supporting Information), and experimental validation using synthetic input combinations of given DPYD genotypes demonstrates expected fluorescence patterns programmed for each recommended dosage (Figure ). Another two linear models were also successfully trained for interpreting SLCO1B1 and UGT1A1 genetic profiles into dose recommendations (Figure , Supporting Information), suggesting that the DNA decoding circuit is highly programmable to guide personalized healthcare. In addition to guiding drug regimes, another rapidly emerging precision medicine tool is the polygenic risk score (PRS), which estimates a person's genetic risk of developing a particular disease by cumulating the effects of various risk alleles with relatively small individual effects on incidence.  Though the clinical validity and utility of polygenic risk profiling for risk stratification, screening, and prediction have been evaluated in detail,  it remains challenging to implement PRS in routine clinical and public health practice,  which prompted us to explore the utility of DNA decoding circuit in the molecular calculation of PRS and PRS‐informed disease screening, considering its high programmability and modularized circuit layout. We thus abstracted a PGS scoring system for pancreatic cancer from the Polygenic Score Catalog ( http://www.pgscatalog.org ).  This PGS scoring system includes 22 risk variants with varying effect sizes (Odds Ratio, OR), and PRS could be calculated as a weighted‐sum score of all effect alleles using their effect sizes from GWAS as weights ( Figure ). Notably, the allelic dosage was defined by the allele copy number (that is, “1” for a heterozygote and “2” for a homozygote), whereas each allele has a binary state of 0 or 1 in our previously engineered circuit. To achieve a ternary input mode, the multiplicative coefficient of effect alleles was programmed to double its original lnOR while the coefficient of other alleles was set to have a negative lnOR (Figure ). In this way, the concentration of weighted species is 2lnOR if the allele is homozygous but reduces to lnOR if the allele is heterozygous, and there are none of positively weighted species in the absence of the allele. Weighted species of each risk locus were then integrated by a summation gate, thereby setting the stage for the molecule calculation of PRS (Figure ). Using synthetic input patterns to simulate homozygous/heterozygous allele samples, we experimentally demonstrated that the DNA circuit could implement reliable molecular calculation of risk score for individual risk locus (Figure , Supporting Information). We also simulated the circuit behavior with the genotypic dataset of 2504 individuals from the 1000 Genomes Project Phase 3 and statistically compared output concentrations of weighted‐sum species with their supposed PRS(s), which were observed to be normally distributed (mean: 4.30; standard deviation: 0.69) (Figure ). Risk stratification based on supposed PRS(s) of the 2504 individuals resulted in three categories: low (bottom decile), intermediate (2nd–8th decile), and high (top decile) (Figure ), and we noted that the stoichiometric quantity of the weighted‐sum output was near to the supposed PRS for most test cases and demonstrated a similar distribution with that of supposed PRS (mean: 4.27; standard deviation: 0.34) (Figure ). Accordingly, we reasoned that high‐risk individual with top10% PRS could be efficiently screened by fluorescent outputs when a comparator module was incorporated into the circuit (Figure ). Using synthetic input patterns of individuals with varying PRS levels (Table , Supporting Information), circuit behaviors were expectedly demonstrated to have an ON trajectory in test cases with top1% PRS and top5% PRS, while maintaining an OFF state for test cases with intermediate PRS (Figure ), suggesting that the DNA‐powered genetic decoder has great potential to efficiently and straightforwardly implement PRS‐informed individualized disease risk prediction and perform populational risk stratification. Compared with the well‐established Sanger sequencing for genetic testing and interpretation (Table , Supporting Information), our DNA‐computation‐based platform shows advantageous features including high multiplexity, automatability, and target expandability, as well as fast turnaround time.   Discussion Herein we present a programmable DNA decoding circuit for straightforward genotype‐phenotype interpretation. We have validated the effectiveness of the DNA‐powered genetic decoder in human genomic DNA samples and demonstrated that the intelligent genetic decoder can be readily reconfigured for diverse genetic services, like individualized medication and cancer risk assessment. To enable parallel molecular computation over a given set of human alleles, a signal‐transforming approach based on thermostable ligase was developed. This approach provides a flexible interface for molecular computing platforms to process biological information by circumventing a series of handing procedures to acquire sequence information of each genetic locus, thereby setting the stage for molecular computation in automated and time‐efficient genetic interpretations. We experimentally demonstrated that this approach could efficiently convert up to five CYP2C19 alleles to recognizable DNA inputs for the circuit in a single reaction. The striking contrast between the thermodynamic favorability of linked and unlinked probes only allows the linked probe to initiate a computing module. Additionally, increasing ligase and primer concentration, as well as extending the annealing time and thermal cycle number, can improve the ligation efficiency. It should be noted that the ligase‐based approach, though validated in the context of single nucleotide variations, is also amenable to transforming other types of genetic variations, such as insertions and deletions, which collectively constitute an overwhelming proportion of human genetic variations.  Molecular programming enables the DNA decoding circuit to implement efficient and straightforward genetic interpretation at the molecular level based on an interpretable linear model, wherein the weights of each allele are sized by respective effects on phenotypic outcomes. That is, a risk allele may obtain a negative weight while a protective allele gets a positive weight. Accordingly, the DNA decoding circuit can make a reliable prediction based on the overall genetic effect and can be readily reprogrammed to satisfy a wide variety of smart applications by varying allelic weight array and/or rewiring circuit modules. Given the modularity of the circuit layout, weight adjustment at the molecular level can be flexibly achieved by simply tunning concentrations of arithmetic gates rather than rewiring the interlinked network of the entire circuit. Also, modularized circuit design enables successful molecular computation of polygenic risk scores, which is rapidly evolving in the academic field but rarely incorporated into healthcare practice, and when combined with the thresholding module, the screening of risky individuals was successfully realized. Further, rationally designed threshold modules help to satisfy multi‐categorical classification requirements in clinical practice, rather than report an “either‐or” answer (i.e., YES/NO) about phenotypical classes and medical decisions with diverse alternatives. Notably, the massively parallel computing ability of the DNA decoding circuit renders it highly scalable with the number of input alleles for more sophisticated tasks, while not markedly increasing the turn‐around‐time. Instead, the workload required for conventional genotyping techniques increases remarkably with the number of genetic variations. For large‐scale DNA strand‐displacement circuits, of particular note is to minimize unspecific crosstalk and circuit leakage, which may be effectively achieved by carefully checking the reactive orthogonality of DNA molecules, introducing a clamp domain into output strands, and purifying computing components. With the high‐density storage and massively parallel computing ability, DNA‐based molecular computation has long been considered a powerful substitute for traditional silicon chips, nevertheless, the most competitive advantage of DNA computation in the healthcare field is its bio‐compatibility, making it a good fit for genetic interpretation at the molecular level. We leveraged this strength and successfully realized the full integration of both genetic detection and data interpretation into a single tube, thereby eliminating additional in silico computation over genetic data acquired by specialized equipment. Though the entire workflow used in this study relies on a PCR machine for pre‐amplification and one‐pot LDR and a microplate reader for the fluorescence readout, however, isothermal amplification methods, for example, recombinase polymerase amplification (RPA), utilize engineered recombinase enzyme instead of thermal denaturing to help primers invade into double‐stranded DNA,  therefore stand as a promising alternative for the isothermal implementation of pre‐amplification and LDR. While fluorescence spectroscopy is conventionally used for the readout of DNA computing systems, there are many emerging alternatives for point‐of‐care readout,  for example, a DNA‐based motor can convert computational results into mechanical motions that can be readily visualized by smartphone camera,  thereby largely reducing the instrument required for the readout of molecular computing platforms. Also of note is that the turn‐around time for the DNA computing system can be dramatically sped up from hours to minutes by a cyclic freeze/thaw approach.  These open up the possibility of the molecular computing system for wide‐scale cell‐free scenarios, including DNA methylation‐based in vitro epigenetic diagnostics, noninvasive prenatal screening, and molecular subtyping for precise therapeutic implications. Further, all step‐by‐step procedures can be integrated into a completely automated workflow assisted by the robotic liquid handling system, revealing the supremacy of our DNA computation platform in high accuracy and automation level.  Given all that, the DNA‐based bio‐computing platform stands as a real‐world demonstration of the utility of molecular computing systems in the intelligent genotype‐phenotype translation of human genomic samples. It can be envisioned to facilitate personalized medication, genetic counseling, and disease risk prediction in primary care and will inspire more powerful and intelligent healthcare applications for improving health outcomes for the general public.   Experimental Section Human Genomic DNA and Synthetic Reference Samples For the preparation of synthetic DNA samples, a 4.81 kb DNA fragment spanning the partial wild‐type CYP2C19 gene was chemically synthesized and cloned into the commercial pUC57 vector, resulting in a plasmid referred to as pUC57‐CYP2C19‐WT. The resulting plasmid then served as a template for generating single nucleotide variants via site‐directed mutagenesis, resulting in another three plasmids carrying mutation at different locations: pUC57‐CYP2C19*2, pUC57‐CYP2C19*3, and pUC57‐CYP2C19*17. After sequence verification, these successfully constructed plasmids were respectively transformed into and propagated in the Escherichia coli Top10 strain to produce a large number of plasmids used as the DNA template in ligase‐dependent reactions. For the preparation of human genomic DNA samples, 500 µL human peripheral intravenous blood samples were used to extract genomic DNA using a blood genomic DNA extraction kit (Tiangen Biotech, DP348). This study was approved by the Ethics Committee at Renji Hospital, School of Medicine, Shanghai Jiao Tong University (RA‐2021‐103). Informed consent forms were obtained from all human participants recruited in this study. CYP2C19 genetic polymorphisms of each sample were verified using Sanger sequencing. Preamplification and SNP Genotyping of Genomic Samples The PCR for the preamplification and SNP genotyping of CYP2C19 with Sanger sequencing was performed using 0.02 units µL −1 Q5 polymerase (New England Biolabs, M0491S), 1× Q5 reaction buffer, 200 n m forward and reverse primers (detailed sequences were listed in (Table , Supporting Information), and 200 µ m dNTPs. The thermal cycling program consisted of 98 °C for 30 s, followed by 30 cycles of 98 °C for 10 s, 65 °C for 30 s, and 72 °C for 20 s. Primer Design for Ligase‐Dependent Reaction Ligase‐dependent reaction (LDR) is used to transform allelic information into corresponding ssDNA inputs for DNA multiplication gates. Each allele has a specific upstream primer which has a particular toehold domain on its 5’ end, and a downstream primer whose 3’ end has been phosphorylated. The upstream allele‐specific primer was expected to be linked with the downstream primer by 9°N ligase. Of note is that the downstream primer is shared by the two alternative alleles at the same locus. To design primer pairs for multiplex alleles, sequences flanking the polymorphic locus of each allele were obtained from NCBI's dbSNP database ( https://www.ncbi.nlm.nih.gov/snp/ ). The primer melting temperature ( T m ) of each primer was then estimated by the online reaction temperature calculator ( https://ligasecalc.neb.com/#!/ligation ). T m of all primers was set to 55 °C, resulting in a list of primer candidates. Next, all candidates of linked primers were subjected to a check of secondary structure at 37 °C by NUPACK.  For linked primers with hairpin structures, intentional mismatches were introduced into the stem domain to make the hairpin structure thermodynamically unfavorable. Note that the base at the 3’ end of all allele‐specific primers should be exactly complementary to the targeted allele. Then, the Tm of each primer candidate with intentionally induced mismatches at the stem domain was then estimated by the reaction temperature calculator ( https://ligasecalc.neb.com/#!/ligation ). Since a single mismatch can influence Tm dramatically, concentrations of primes with a single mismatch were increased to ensure the hybridization efficiency of the primers to the target template. The incubation temperature of all primer pairs should be within a similar range (50 °C ± 3 °C). Ligase‐Dependent Reaction Ligase‐dependent reaction were performed using 1.6 units µL −1 9°N ligase (New England Biolabs, M0238S), 1 × 9°N ligase buffer, 5 n m plasmid DNA and all necessary primer pairs. Detailed sequences of all primers were specified in Table , Supporting Information. The thermal cycling program consisted of 94 °C for 5 min, followed by 30 cycles of 94 °C for 10 s and 50 °C for 90 s. The volume of LDR product used for 100 µL DNA circuit was 20 µL. Sequence Design For all sequences, no more than four A or T and no more than three C were used in a row. Use A, C, and T only for all output strands to reduce secondary structures and undesired interactions; Keep GC content at 30–70% for gate strands. Any two sequences in a sequence pool should not share more than 35% of the domain length. Purification of DNA Duplexes DNA oligonucleotides were chemically synthesized and purified by Sangon Biotech. All unmodified DNA oligonucleotides were purified using ULTRAPAGE, and modified DNA oligonucleotides were purified using HPLC. All oligonucleotides were suspended at 100 µ m in 1× TE buffer with 12.5 m m MgCl 2 and stored at −20 °C. Top and bottom strands of DNA duplexes were pre‐mixed in a 1:1 ratio and annealed at 20 µ m in a thermal cycler by heating to 95 °C for 5 min and then ramping down to 25 °C at a rate of 1 °C per 1 min. Purification of DNA duplexes was subsequently performed using 12% native polyacrylamide gel electrophoresis, and then double‐stranded bands were carefully cut out and transferred into 1× TE buffer with 12.5 m m MgCl 2 to allow free diffusion of DNA out of the gel. Absorbance at 260 nm was determined by NanoDrop One (Thermo Fisher) for quantitation of DNA duplex concentration using extinction coefficients. All purified DNA gates were stored at 4 °C for later use. Detailed sequences of DNA strands used in the circuit were listed in Table , Supporting Information. Fluorescence Kinetics Kinetic fluorescence measurements were performed using the Synergy H1 Hybrid Multi‐Mode Reader (BioTek). Unless otherwise specified, all reactions were carried out in 1× TE buffer with 12.5 m m Mg 2+ and kept at 37 °C throughout the reaction. Excitation/Emission wavelengths were set to 493/522 nm for FAM, 535/575 nm for VIC, 585/615 nm for ROX, and 643/671 nm for Cy5. The master mixture containing all circuit components, except input strand(s) or LDR products, was initially added to the 96‐well all‐black plate to record the background fluorescence, followed by the addition of input strand(s) or LDR products. Fluorescence signals were thereafter collected every 2 min for each well. Statistical Analysis Raw fluorescence signals were normalized into the standard concentration of respective output strands by performing linear fitting based on a calibration curve constructed using the steady‐state fluorescence levels of reaction mixtures with a series of reporter initiator strands of predetermined concentrations as input strand(s). For normalization of raw fluorescence intensities of a given fluorophore into relative fluorescence unit (R.F.U), the baseline fluorescence obtained from the initial measurement of the negative control (no input strands) was treated as the minimum (relative fluorescence unit = 0), while the highest fluorescence level in each set of parallel kinetic fluorescence experiments or fluorescence from the standard reporter initiating strands was treated as the maximum (relative fluorescence unit = 1). Statistical analysis of PRS distribution was performed using RStudio (2022.02.0).   For the preparation of synthetic DNA samples, a 4.81 kb DNA fragment spanning the partial wild‐type CYP2C19 gene was chemically synthesized and cloned into the commercial pUC57 vector, resulting in a plasmid referred to as pUC57‐CYP2C19‐WT. The resulting plasmid then served as a template for generating single nucleotide variants via site‐directed mutagenesis, resulting in another three plasmids carrying mutation at different locations: pUC57‐CYP2C19*2, pUC57‐CYP2C19*3, and pUC57‐CYP2C19*17. After sequence verification, these successfully constructed plasmids were respectively transformed into and propagated in the Escherichia coli Top10 strain to produce a large number of plasmids used as the DNA template in ligase‐dependent reactions. For the preparation of human genomic DNA samples, 500 µL human peripheral intravenous blood samples were used to extract genomic DNA using a blood genomic DNA extraction kit (Tiangen Biotech, DP348). This study was approved by the Ethics Committee at Renji Hospital, School of Medicine, Shanghai Jiao Tong University (RA‐2021‐103). Informed consent forms were obtained from all human participants recruited in this study. CYP2C19 genetic polymorphisms of each sample were verified using Sanger sequencing.   The PCR for the preamplification and SNP genotyping of CYP2C19 with Sanger sequencing was performed using 0.02 units µL −1 Q5 polymerase (New England Biolabs, M0491S), 1× Q5 reaction buffer, 200 n m forward and reverse primers (detailed sequences were listed in (Table , Supporting Information), and 200 µ m dNTPs. The thermal cycling program consisted of 98 °C for 30 s, followed by 30 cycles of 98 °C for 10 s, 65 °C for 30 s, and 72 °C for 20 s.   Ligase‐dependent reaction (LDR) is used to transform allelic information into corresponding ssDNA inputs for DNA multiplication gates. Each allele has a specific upstream primer which has a particular toehold domain on its 5’ end, and a downstream primer whose 3’ end has been phosphorylated. The upstream allele‐specific primer was expected to be linked with the downstream primer by 9°N ligase. Of note is that the downstream primer is shared by the two alternative alleles at the same locus. To design primer pairs for multiplex alleles, sequences flanking the polymorphic locus of each allele were obtained from NCBI's dbSNP database ( https://www.ncbi.nlm.nih.gov/snp/ ). The primer melting temperature ( T m ) of each primer was then estimated by the online reaction temperature calculator ( https://ligasecalc.neb.com/#!/ligation ). T m of all primers was set to 55 °C, resulting in a list of primer candidates. Next, all candidates of linked primers were subjected to a check of secondary structure at 37 °C by NUPACK.  For linked primers with hairpin structures, intentional mismatches were introduced into the stem domain to make the hairpin structure thermodynamically unfavorable. Note that the base at the 3’ end of all allele‐specific primers should be exactly complementary to the targeted allele. Then, the Tm of each primer candidate with intentionally induced mismatches at the stem domain was then estimated by the reaction temperature calculator ( https://ligasecalc.neb.com/#!/ligation ). Since a single mismatch can influence Tm dramatically, concentrations of primes with a single mismatch were increased to ensure the hybridization efficiency of the primers to the target template. The incubation temperature of all primer pairs should be within a similar range (50 °C ± 3 °C).   Ligase‐dependent reaction were performed using 1.6 units µL −1 9°N ligase (New England Biolabs, M0238S), 1 × 9°N ligase buffer, 5 n m plasmid DNA and all necessary primer pairs. Detailed sequences of all primers were specified in Table , Supporting Information. The thermal cycling program consisted of 94 °C for 5 min, followed by 30 cycles of 94 °C for 10 s and 50 °C for 90 s. The volume of LDR product used for 100 µL DNA circuit was 20 µL.   For all sequences, no more than four A or T and no more than three C were used in a row. Use A, C, and T only for all output strands to reduce secondary structures and undesired interactions; Keep GC content at 30–70% for gate strands. Any two sequences in a sequence pool should not share more than 35% of the domain length.   DNA oligonucleotides were chemically synthesized and purified by Sangon Biotech. All unmodified DNA oligonucleotides were purified using ULTRAPAGE, and modified DNA oligonucleotides were purified using HPLC. All oligonucleotides were suspended at 100 µ m in 1× TE buffer with 12.5 m m MgCl 2 and stored at −20 °C. Top and bottom strands of DNA duplexes were pre‐mixed in a 1:1 ratio and annealed at 20 µ m in a thermal cycler by heating to 95 °C for 5 min and then ramping down to 25 °C at a rate of 1 °C per 1 min. Purification of DNA duplexes was subsequently performed using 12% native polyacrylamide gel electrophoresis, and then double‐stranded bands were carefully cut out and transferred into 1× TE buffer with 12.5 m m MgCl 2 to allow free diffusion of DNA out of the gel. Absorbance at 260 nm was determined by NanoDrop One (Thermo Fisher) for quantitation of DNA duplex concentration using extinction coefficients. All purified DNA gates were stored at 4 °C for later use. Detailed sequences of DNA strands used in the circuit were listed in Table , Supporting Information.   Kinetic fluorescence measurements were performed using the Synergy H1 Hybrid Multi‐Mode Reader (BioTek). Unless otherwise specified, all reactions were carried out in 1× TE buffer with 12.5 m m Mg 2+ and kept at 37 °C throughout the reaction. Excitation/Emission wavelengths were set to 493/522 nm for FAM, 535/575 nm for VIC, 585/615 nm for ROX, and 643/671 nm for Cy5. The master mixture containing all circuit components, except input strand(s) or LDR products, was initially added to the 96‐well all‐black plate to record the background fluorescence, followed by the addition of input strand(s) or LDR products. Fluorescence signals were thereafter collected every 2 min for each well.   Raw fluorescence signals were normalized into the standard concentration of respective output strands by performing linear fitting based on a calibration curve constructed using the steady‐state fluorescence levels of reaction mixtures with a series of reporter initiator strands of predetermined concentrations as input strand(s). For normalization of raw fluorescence intensities of a given fluorophore into relative fluorescence unit (R.F.U), the baseline fluorescence obtained from the initial measurement of the negative control (no input strands) was treated as the minimum (relative fluorescence unit = 0), while the highest fluorescence level in each set of parallel kinetic fluorescence experiments or fluorescence from the standard reporter initiating strands was treated as the maximum (relative fluorescence unit = 1). Statistical analysis of PRS distribution was performed using RStudio (2022.02.0).   The authors declare no conflict of interest.   Supporting Information Click here for additional data file.   