Stepwise In Silico Specificity Approach for NGT PCR-Based Detection Methods in the Context of the EU GMO Regulations

1. Introduction

Genetically modified organisms (GMOs) are subject to regulations in Europe and worldwide through domestic legislation and the Cartagena Protocol on Biosafety (CPB), an international GMO treaty by the United Nations. In Europe, prior to market approval of GMOs and derived food and feed products, GMO event-specific methods for their detection, identification, and quantification need to be in place according to European Commission (EC) Regulation No 1829/2003 on genetically modified food and feed. The analytical framework to detect GMOs in Europe relies on validated qualitative and quantitative real-time PCR methods, forming part of the single harmonized, time-limited, and transparent procedure for the authorization of GM food and feed in the European Union (EU) (EC,2021).

However, recent developments in the field of genetic engineering, known as ‘new genomic techniques’ (NGT), have raised concerns over the suitability of such detection methods and the possibility of implementing current GMO regulations due to small genomic alterations are an umbrella term used to describe various techniques that can alter the genetic material of an organism and have emerged or been developed since 2001 when the existing genetically modified organism (GMO) legislation was adopted (EC, 2021). These techniques are based on site-directed nucleases (SDNs) (GROHMANN et al., 2019) and encompass different enzymes, such as meganucleases (LAGLIDADG endonucleases; EMN), zinc finger nucleases (ZNF), effector transcription factor nucleases (TALENs), and Clustered Regularly Interspaced Short Palindromic Repeats system (CRISPR-Cas9) (SPRINK et al., 2016) (BROLL; BRAEUNING; LAMPEN, 2019) (METJE-SPRINK et al., 2019) (GAJ; GERSBACH; BARBAS, 2013) (OSAKABE; OSAKABE, 2015).

In 2018, the Court of Justice of the European Union (CJEU) issued a significant decision stating that gene-editing techniques, including CRISPR-Cas9, should be classified and regulated as GMOs. As a result, any techniques that modify the genetic material of plants or animals using NGT should be subject to the same regulations applied to GMOs. In this context, the EU provide significant steps to implement legislation ensuring the traceability of GMO-containing foods and derivatives in the market, providing scientific guidelines and technical documents to empower consumers with informed choices through proper labeling (ENGL, 2008) (EUROPEAN LABORATORIES, 2015) (ARAKI; NOJIMA; ISHII, 2014) (FRAITURE et al., 2015).

In the context of detecting NGTs, two main steps are crucial for the development of PCR-based analysis methods: i) in-silico tests and ii) empirical tests. The in-silico computational assessments involve the verification of the compatibility of primer sequences and avoid unintended amplifications of other genomic regions or sequences from different organisms. Ensuring the specificity of primers designed for the edited organism is essential, taking into consideration any base differences. PCR relies on the complementarity between primers and the target DNA for successful amplification. The design of specific primers and probes involves predicting secondary structures to avoid interferences, such as primer dimers or hairpin formations, which can affect amplification specificity. Additionally, evaluating sequence similarity with other organisms ensures that the target sequences of interest do not share similarities with sequences from different organisms (ENGL, 2021). The technical complexity of designing a PCR system capable of distinguishing specific sequences with accuracy and reliability requires both in-silico and experimental testing. Some studies have successfully distinguished NGT canola and rice varieties using PCR (WEIDNER et al., 2022) (FRAITURE et al., 2022) (ZHANG et al., 2022).

To address the challenges associated with detecting edited organisms, our approach can be applied to explore and characterize NGT-modified organisms, as it adheres to the principles of open science, ensuring transparency in data and efforts for the detection, identification, and biosafety of these organisms. In this study, we present a step-by-step bioinformatics analysis to predict primer specificity requirements in the context of the EU GMO law and regulations, employing two NGT model sequences from Arabidopsis thaliana, each harboring distinct mutation patterns (single polymorphism for genotype 1 and a 3 bp deletion for genotype 2) in the target growth-regulating factor (GRF) gene. We expect our approach to offer a theoretical foundation for assessing primer specificity, complemented by empirical tests to support the hypothesis of selective target amplification.

2. Materials and Methods

2.1. Overall Description of the Stepwise Approach

The theoretical basis of real-time polymerase chain reaction (RT-PCR) was used to screen mutations delivered by a CRISPR-Cas9 technique at the target site. The first step identified in our approach is the sanger sequencing of mutated alleles for confirmatory purposes. After confirmation of mutation pattern in selected genotypes, the next step is a comprehensive search on existing databases for sequence similarities and the natural variability around the mutated sequence. The outcomes of the search analysis were manually evaluated for the presence of mismatches in regions corresponding to the amplicon sequence. Lastly, an in-silico PCR testing is performed to estimate the amplification rates for each PCR system and a second confirmatory steps for primer and probe specificity. Each of these steps is described in detail in the next sections and a full picture of the overall approach is presented below (Figure 1). We performed all the steps in a model case for two gene-edited genotypes of Arabidopsis. The description of the biological material used in this study is provided in the next section.

2.2. Confirmatory Sequencing Analysis of CRISPR/Cas9 Mutations

The seeds from two CRISPR-Cas9 genotypes, grf1-3 and grf8-61, as well as the control of Arabidopsis thaliana, were obtained from the Eurasian Arabidopsis Stock Center (NASC). The transformation of plants was carried out using the Agrobacterium strain GV3101 through the floral dip protocol. In this method, Arabidopsis thaliana with a Columbia background served as the wild-type strain for transformation and CRISPR/Cas9 mutagenesis. The plants were immersed in an Agrobacterium suspension to facilitate the transfer of CRISPR/Cas9 components into the plant cells (ANGULO et al., 2020). The grf1-3 genotype is a null mutant and transgene-free. The mutation is located in the GRF1 (Growth-Regulating Factor 1) gene, specifically in the locus AT2G22840. This mutation involves a single nucleotide polymorphism (SNP) with a Guanine insertion at position 10048382 on the positive strand of chromosome 2. On the other hand, the grf8-61 genotype carries a three-base pair (bp) deletion, resulting in the mutation of the GRF8 (Growth-Regulating Factor 8) gene located in the locus 2134986.Additionally, the Columbia ecotype, sequenced in the Arabidopsis Genome Initiative (https://www.ebi.ac.uk/ols/ontologies/efo/terms?short_form=EFO_0005147), was used in the study. Three individual plants per genotype were germinated in sterilized soil placed in pots. After 10 days of incubation at 4°C in the dark, the plants were transferred to a greenhouse with a light period of 8 hours/day and an intensity of approximately μmol m-2 s-1 at 25°C for growth. At the growth stage (G), ten rosette leaves per plant were used for DNA isolation using the Plant DNAzolTM Reagent from Invitrogen Co., US. The target regions from genotypes grf1-3, grf8-61, and Columbia were amplified through conventional PCR using Sanger primers with the following cycling conditions: 2 minutes at 94°C, followed by 35 cycles of 45 seconds at 94°C, 30 seconds at 60°C, 1 minute at 72°C, and a final extension of 10 minutes at 72°C (Applied Biosystems™ Veriti™ Thermal Cycler). The PCR products from each genotype were purified and sequenced using the Applied Biosystems Genetic Analyzer (3500×, São Paulo, Brazil). The obtained Sanger sequences from each genotype were checked for quality using Chromas version 2.6.6 (Technelysium Pty Ltd., Brisbane, Australia). The alignment comparison with default settings (MUSCLE) was performed to compare the sequences obtained from each Sanger fragment to the GRF1 (AT2G22840; ID: 816815; Chr2) and GRF8 (AT4G24150; ID: 828515; Chr4) genes of Arabidopsis thaliana (NCBI GenBank taxid:3702).

2.3. Amplicon Sequence Search for Natural Variants

The grf1-3 and grf8-61 amplicon sequences and their wild type of counterpart, Columbia (Col-01), were searched against the National Center for Biotechnology Information database (NCBI) using the Basic Local Alignment Search Tool (BLAST). The databases accessed were nucleotides (nr/nt) available at GenBank + EMBL + DDBJ + PDB + RefSeq databases. Patent strings, phase strings 0, 1 and 2 HTGS, EST, STS, GSS, WGS, TSA and greater than 100MB were excluded. To be considered as valid hits, outcome sequences retrieved through BLAST needed to follow the 80% coverage parameter; 80 to 98% identity level, and sequence hits below 10 bp of mismatches were considered eligible. The presence of mismatches in regions corresponding to the primer and probe sets were manually verified using the Multiple Sequence Alignment Viewer (MSA) graphic displays for NCBI nucleotide alignment. The position of the single point mutation is indicated in Figure 3A) and the positions of the exclusions, indicated Figure 2B) were used to search for position polymorphisms in a database composed of 1,001 variant genomes of Arabidopsis thaliana (https://1001genomes.org/). The reason behind the strategy to use amplicon, instead of classical primer and probe sets, was to avoid false positive hits when only one primer set, or probe matches the sequence, but no effective amplification is expected at similar amplicon size.

2.4. RT-PCR Primer and Probe Design

The design of primers and probes for detecting CRISPR-Cas9 mutations in the grf1-3 and grf8-61 genotypes of Arabidopsis thaliana was performed using the Primer-3 Plus software (ROZEN; SKALETSKY, 2000) (UNTERGASSER et al., 2007). To ensure primer specificity and avoid PCR competition, the sequences of each primer set were manually verified in both DNA strand orientations to contain the CRISPR-Cas9 mutation. First, the size of the amplicon for the single-copy reference gene was taken as a basis to design the amplicons for each edited event (Figure 2). Subsequently, the selection of primer and probe sets was based on specific parameters. All oligonucleotides were checked for the GC content of the primers (between 45-55%) and the melting temperature (Tm), which was set 3 to 5°C below the Tm of the probe (between 62 and 68°C) (ROBLEDO; BEGGS; BENDER, 2005). Additionally, secondary structures and hairpins were considered (MIQE, 2009) (ENGL, 2021). The selected primers were also evaluated for similarity to the endogenous reference gene amplicon (AT1G03400) of Arabidopsis thaliana (YANG et al., 2018).

2.5. In-Silico PCR Testing

In addition, a second confirmatory step was performed for in-silico PCR testing using PrimerBLAST software (Ye et al., 2013). This analysis aims to screen for primer annealing in directions from 5’ to 3’ between all combinations, including forward-to-reverse, forward-to-forward and reverse-reverse. The in-silico PCR testing was performed in all available species in the database and no amplification was detected in data set species. For each amplicon template (grf1-3 and grf8-61) the PCR products were limited to the size of 100 to 300 bp, and the average melting temperature (Tm) (Supplementary File S2) was adjusted. The database analyzed was from genomes Eukariotic For the stringency of the primers, it was considered that the primer should contain at least 2 unintended targets, including at least 2 incompatibilities in the last 5 bps at the 3’ end. The number of target sequences retrieved for each template was considered an indicative to the reproducibility and chances of finding amplification from different target sequences than the intended one.

3. Results

3.1. Sequence Confirmation of Mutated Alleles

After sanger sequencing of mutated genotypes, and wild type non-mutated genotype, intend mutations were confirmed for each new allele obtained by CRISPR-Cas9 (Figure 2).

Figure 2. Alignment between the sequences obtained through the sequencing data. Panel A shows the comparison between the sequences obtained from the Columbia (Col-1) and the grf1-3 CRISPR-Cas9 allele. The location of the Guanine (G) insertion is demonstrated in the respective position of Arabidopsis thaliana genome. Panel B shows the comparison between sequences (Col-1) and the gfr8-61 CRISPR-Cas9 allele. The three bases deleted, Guanine (G) Cytosine (C) and Adenosine (A) are shown in the same respective genome. For both cases, the site at 3 bp after the PAM Cas9 was useful for mutations localization.

Figure 2. Alignment between the sequences obtained through the sequencing data. Panel A shows the comparison between the sequences obtained from the Columbia (Col-1) and the grf1-3 CRISPR-Cas9 allele. The location of the Guanine (G) insertion is demonstrated in the respective position of Arabidopsis thaliana genome. Panel B shows the comparison between sequences (Col-1) and the gfr8-61 CRISPR-Cas9 allele. The three bases deleted, Guanine (G) Cytosine (C) and Adenosine (A) are shown in the same respective genome. For both cases, the site at 3 bp after the PAM Cas9 was useful for mutations localization.

3.3. In-Silico Specificity Assessment

The search for amplicon sequences for wild type GRF1 gene in the NCBI nucleotide collection (nr/nt) recovered a total of 120 hits/occurrences in 39 organisms. Of these, 78 hits and 15 organisms belong to the Brassicales Order (Supplementary Archive S3). We recovered 58 hits related to the genus Arabidopsis, of these, 56 accesses related to Arabidopsis thaliana (100% identity) and the other 2 hits belong to the species-genus Arabidopsis lirata (99.23% identity) and Arabidopsis arenosa (98.46% identity). In five hits, the sequences were 100% identical to the GRL1 gene of the growth-regulating protein (B2CU94_ARATH). About the wild type GRF8 sequence, 57 hits were recovered with 100% coverage and 100% identity in Arabidopsis thaliana (Supplementary File S4). Of these, 25 accesses were 100% identical to the GRL8 gene of the protein growth regulatory factor in Arabidopsis thaliana (B2CUI8_ARATH). These results show that no significant homology is found for the two gene alleles outside the Arabidopsis genus (Supplementary File S5).

For amplicon grf1-3 CRISPR-Cas9 (mutated), 36 hits with sequences with 100% coverage were found. Two Arabidopsis species, A. Lirata and A. Arenosa, presented an identity of 99.71% and 98.47%, respectively to the mutated GRF1amplicon. Both species showed identical sequence with a SNP variation in position 15. Additionally, A. arenosa presents a mismatch at position 108. Other species, such as Camelina sativa and Camelina hispida showed 100% coverage but lower identity levels reaching 98.48% and 97.71%, respectively due to different mismatches. We also found 10 hits from 4 species from Brassica genus which presented 95.5% identity. Raphanus sativus species also showed 2 sequence hits with coverage and identity of 98% and 96.12%, respectively. No other polymorphisms were found in the database for Arabidopsis thaliana for the given amplicons. Matching results are presented in Table 1 for grf1-3.

For the grf8-1 amplicon, 26 hits were recovered with 96% coverage due to the presence of gaps in the initial bases, corresponding to (C, A, G, C, T) the forward primer. After the 6-mer position, the sequences showed 100% coverage. In two of those hits, the identity was 97.54% due to the presence of mismatches in the probe sequence. A 99.20% identity and 98% coverage was also found for two other hits that have a mismatch in the 5’ 5-mer a transversion-type mutation in DNA (A) position 12538392 in the forward primer. No other sequences were found with 100% alignment for this amplicon. The variation across these hits indicates that this 5pb can be composed of polymorphic nucleotides. However, the three excluded bases were not found in the variant database for Arabidopsis thaliana. Table 3 complies the results obtained for the grf8-61 amplicon search against the database.

3.4. In-Silico PCR Performance

A pair of primers should only amplify the target sequence. This is especially challenging when the target region differs by one or a few nucleotides from other potential targets. The in-silico PCR results presented 57 potential hits for the grf1-3 PCR system. Of this, 39 hits were observed for the forward primer containing a Cytosine (C) base inserted as a mismatch. The other 18 hits were distributed in Arabidopsis hispida (1 hit) and Arabidopsis arenosa (1) and Brassicas sp (16). Despite having mismatches in the amplicon against Brassica napus, Brassica oleracea, Brassica Rapa, these species were still considered a potential amplification sequence in our in-silico PCR prediction. These results suggest that the primer set designed to distinguish grf1-3 CRISPR-Cas9 genotype has the potential to amplify sequences in other species, such as Arabidopsis and Brassica.

Table 3. Summary of in silico-PCR amplification for the grf1-3 amplicon PCR in silico matrix for the grf1-3 amplicon. These theoric condition indicate unintended targets. Thus, primers may not be specific for the PCR template.

Total of mismatches	Number/hits analyzed	Sequences corresponding perfectly to the primer	Number of Blast hits recovered	Possible discrimination between the grf1-3 genotype and other lines
1	39	0	57	Low
2	17	0
3	1	0

Table 3. Summary of in silico-PCR amplification for the grf1-3 amplicon PCR in silico matrix for the grf1-3 amplicon. These theoric condition indicate unintended targets. Thus, primers may not be specific for the PCR template.

Total of mismatches	Number/hits analyzed	Sequences corresponding perfectly to the primer	Number of Blast hits recovered	Possible discrimination between the grf1-3 genotype and other lines
1	39	0	57	Low
2	17	0
3	1	0

The primer set design for the detection of 3bp deletion in the Growth Factor Regulator 8 gene showed a unique target in Arabidopsis thaliana. Therefore, the primer set designed can be considered specific for the grf8-61 genotype after consideration of the in silico PCR prediction.

Table 4. Summary of in-silico PCR amplification result for the grf8-61. The amplicon correctly identifies the new CRISPR-Cas9 grf8-61 genotype. The primers demonstrate specificity in-silico PCR.

Total of mismatches	Number/hits analyzed	Sequences corresponding perfectly to the primers	Number of Blast hits recovered	Possible discrimination between the grf8-3 genotype and other lines
1	1	1	1	High

Table 4. Summary of in-silico PCR amplification result for the grf8-61. The amplicon correctly identifies the new CRISPR-Cas9 grf8-61 genotype. The primers demonstrate specificity in-silico PCR.

Total of mismatches	Number/hits analyzed	Sequences corresponding perfectly to the primers	Number of Blast hits recovered	Possible discrimination between the grf8-3 genotype and other lines
1	1	1	1	High

4. Discussion

DNA breaks caused by site-directed nuclease, such as CRISPR-Cas, in genetically modified organisms (GMOs) can be repaired through two major pathways, one pathway joins the ends from a template (HR) and another pathway that joins the ends without homology (NHEJ). In plants produced by NBTs, predominantly, the NHEJ pathway is used to obtain knockout of genes and mutants with characteristics of interest to agriculture (YAN et al., 2020) (Rodríguez-Leal et al., 2017) (Sedeek et al., 2019).

Mutants generated by this repair mechanism have one or a few InDels other than their unmutated part, and these different nucleotides are used as identifiers to characterize the new GMO. In this study, two transcriptional growth regulatory genes (GRFs) regions in Arabidopsis thaliana were used as case studies to develop a stepwise approach for in silico prediction of specific PCR primer sets. The GRF protein gene family is involved in growth and development, and the stress response (Liu et al., 2016) (Huang et al., 2021, Li et al., 2021) (OMIDBAKHSHFARD et al., 2015).

The in silico screening of the amplicon sequences of the native genes of A. thaliana Columbia in the NCBI nucleotide database allowed us to identify that the primer sequence used in the GRF1 gene has an identity of 99% to Arabidopsis lyrate, 97.24% to Arabidopsis arenosa, and 94.6% to agricultural species Brassica napus, Brassica rapa, Brassica rapa subsp. Pekinensis and Brassica oleracea var oleracea. Another 15 organisms belonging to the order Brassicales were found with lower similarity, above 6 bp of mismatches. The amplicon sequence search for GRF8 gene showed only one hit in Arabidopsis thaliana with 100% identity, thus showing our own native sequence.

As expected, the BLAST of amplicon grf1-3 indicated that the single nucleotide polymorphism is the only difference in 35 sequence hits with 99.8% similarity. Similar results were found by (WEIDNER et al., 2022) analyzing amplicon of primer sets for the detection of gene-edited canola. In their study, the BLAST search yielded hits in B. oleracea and B. carinata, in which the only difference between the targets was the last nucleotide at position 22 – in the position of the complementary reverse primer to the base for detection (WEIDNER et al., 2022).

Most gene edited organisms will contain mutations in conserved regions and/or relevant exons for trait accomplishment. This is the case for gene edited Camelina sativa and Brassica napus obtained by CRISPR-Cas9 for high oleic acid content (Jiang et al., 2017; Morineau et al., 2017; Oku-Zaki et al., 2018), and the fungus Verticillium longisporum resistant Arabidopsis thaliana and Brassica napus which contain a unique base transition in the CRT1a gene (calreticulin) (PRÖBSTING et al., 2020). From a detection perspective, these GMOs need to show unique amplification products for which the method was designed (EFSA 2017). It should simultaneously distinguish products from naturally occurring variants, or from products obtained by conventional mutagenesis (ENGL, 2019). While real-time PCR amplification is influenced by the primer size, SNP position, hybridization, Taq polymerase fidelity and other PCR conditions; stable non-refracted thermodynamic energy to enhance specificity hybridization can be achieved by the blocking the last 2 pb of 3’ end in strand in primer sequences (HUANG; ARNHEIM; GOODMAN, 1992) (AYYADEVARA; THADEN; SHMOOKLER REIS, 2000).

Whereas grf1-3 and grf8-61 genotypes were exclusive variants, meaning they do not have polymorphic position variants found in the Arabidopsis thaliana 1001 variants database, the results of the in-silico PCR were district for each genotype. The primer set designed for genotype grf1-3 does not distinguish the mutant event from 39 other sequence hits observed in Arabidopsis arenosa, Arabidopsis lyrata, Camelina sativa and Brassica sp. This means, for example, that samples from food mixtures of Arabidopsis thaliana and Camelina sativa, the primer set may not be able to differentiate the modified event of CRISPR-Cas9 for the grf1-3 event. For genotype grf8-61, the results show that the primer set is specific for the detection of the new genotype produced by CRISPR-Cas9 considering the NCBI genetic database search.

Quantitative real-time RT-qPCR method based on hydrolysis probe chemistry is the most widely used methods for GMO detection in Europe (ENGL 2015). In previous studies, amplicon and primer sets were used to identify whether the single point variation in a gene that carries a single nucleotide (adenosine) insertion in rice variety OsMADS26 (locus: Os08g02070). The specificity of the method in distinguishing the modification introduced by gene editing technology using its single variation point was achieved in a method based on 2-plex digital droplet PCR (FRAITURE et al., 2022). Similarly in rice, the identification of genotypes with InDels from 1 to 18 base pairs (bps) for the chlorophyll oxygenase 1 (CAO1) gene region was empirically demonstrated by real-time quantitative (qPCR) and digital droplet (ddPCR) PCR. In both PCR systems, it was possible to estimate genome-edited ingredient content (without in silico approach) in gene edited samples (ZHANG et al., 2021).

However, even when these PCR systems are able to distinguish gene-edited and unmodified parental lines empirically, there is a need to test whether the system is suitable in identifying such GMO events in food mixtures containing several other potential sequence variants. In this case, the value of an in silico approach is to anticipate potential sequence targets already available in genomic databases. A recent controversy shows how relevant it is to perform such prediction analysis. In an earlier study by Chhalliyil and co-authors (2022), empirical data showed specificity of a primer set based on locked nucleic acids (LNA) strategy for detecting and identifying the first commercialized genome-edited plant, the Cibus canola, containing two mutated genes AHAS1C and AHAS3A. Later that same year, Weidner and co-authors (2022) showed that the method might not be specific to the GMO event in case but could also show spurious amplification in other Brassica species.

Whereas the empirical PCR testing is the ultimate proof of method specificity, it is the in silico prediction that is able to anticipate failure to distinguish amplicons from different organisms. Therefore, our study provides a stepwise approach to search for amplicon and primer sets specificity in available databases. Additionally, it allows the elaboration of strategies to increase the amplification specificity in RT-qPCR and/or ddPCR.

Arabidopsis thaliana is a well-sequenced genomic plant model with a broad degree of knowledge of natural variants. However, the search and alignment approach proposed here, as well as in silico PCR testing, identified that other similar sequence occurrences are present at the NCBI-GenBank for the genotype by grf1-3. Therefore, wet lab laboratory testing will have to consider chemical modifications in the PCR system to increase specificity, such as locked nucleic acid technology (LNA), RNAse-H systems, etc (BROCCANELLO et al., 2018) and/or ddPCR startegies to provide specificity in the detection of mutant and wild type variants simultaneously (Fraiture et al., 2022). Finally, considering that most of the genes of interest can be congruent in different species and that in the future NGT plants will be present in food and feed matrices, it is essential to generated information and also systematize and compile information to be included in GMO databases such as Eugenius, and/or JCR to promote the identification of specific sequences relating these GMO events in the future.

5. Conclusions

The findings suggest that primer set designed to identify the grf1-3 genotype containing a single nucleotide polymorphism might generate amplicons in other species. The DNA sequence mismatches that are present are not sufficient to be considered a event-specific PCR system. Under these conditions, the primer pair does not discriminate the grf1-3 genotype of the species of Arabidopsis thaliana itself, other species of the same genus (A. lirata and other phylogenetically close species (i.e., C. Sativa and Brassica sp). Differently, the primer set for grf8-61 genotype is shown to be sufficiently specific for identifying and discriminating the event modified by CRISPR-Cas9.