Effect of endogenous reference genes on digital PCR assessment of genetically engineered canola events.

Droplet digital PCR (ddPCR) has been used for absolute quantification of genetically engineered (GE) events. Absolute quantification of GE events by duplex ddPCR requires the use of appropriate primers and probes for target and reference gene sequences in order to accurately determine the amount of GE materials. Single copy reference genes are generally preferred for absolute quantification of GE events by ddPCR. Study has not been conducted on a comparison of reference genes for absolute quantification of GE canola events by ddPCR. The suitability of four endogenous reference sequences (HMG-I/Y, FatA(A), CruA and Ccf) for absolute quantification of GE canola events by ddPCR was investigated. The effect of DNA extraction methods and DNA quality on the assessment of reference gene copy numbers was also investigated. ddPCR results were affected by the use of single vs. two copy reference genes. The single copy, FatA(A), reference gene was found to be stable and suitable for absolute quantification of GE canola events by ddPCR. For the copy numbers measured, the HMG-I/Y reference gene was less consistent than FatA(A) reference gene. The expected ddPCR values were underestimated when CruA and Ccf (two copy endogenous Cruciferin sequences) were used because of high number of copies. It is important to make an adjustment if two copy reference genes are used for ddPCR in order to obtain accurate results. On the other hand, real-time quantitative PCR results were not affected by the use of single vs. two copy reference genes.

Introduction

Digital PCR is being widely used for the detection and quantification of genetically engineered (GE) events [[1], [2], [3]]. Specific and single copy endogenous reference genes are preferred for absolute quantification of GE events by PCR. Different endogenous reference genes have been used for real-time PCR detection and quantification of GE canola events. Five endogenous reference genes [acetyl-CoA carboxylase (BnACCg8), phosphoenolpyruvate carboxylase (PEP), oleoyl hydrolase (FatA)), high-mobility group protein I/Y (HMG-I/Y) and cruciferin A (CruA)] were compared for specific real-time PCR detection and quantification of Brassica napus [4]. Two different sequences were reported for the cloned fragments of HMG-I/Y, PEP and CruA, indicating the presence of genes in two copies [4]. On the other hand, HMG-I/Y was reported to be a single copy reference gene that can be used for quantification of GE canola events by real-time PCR [5]. It was also reported that the five endogenous reference genes mentioned above were not suitable for real-time PCR quantification of GE canola events as they were not specific between different species and also not stable among cultivars [4]. However, endogenous reference genes such as CruA and HMG-I/Y have been widely used for real-time PCR quantification of GE canola events and acceptable results were reported using the reference genes [examples: [[6], [7], [8]]. Recently, Acyl-ACP thioesterase (FatA(A)) gene was reported to be specific to the A genome of cultivated canola quality oilseed rape (B. napus, B. rapa and B. juncea) and recommended to be used as endogenous reference gene for real-time PCR [9]. There is no published information on the comparison of endogenous reference genes for absolute quantification of GE canola events by droplet digital PCR (ddPCR). The objectives of the study were to: 1. assess suitability of endogenous Cruciferin (Ccf), CruA, FatA(A) and HMG-I/Y reference genes for absolute quantification of GE canola events using ddPCR; 2. determine the influence of DNA extraction methods, DNA quality and cultivar variation on copy numbers of reference genes; and 3. compare the effect of reference genes on ddPCR and real-time PCR results.

Materials and methods

Seed sources

Seeds of Armor BX (OXY235 canola event), Innovator (HCN92 canola event), Legend (non-GE canola), AC Parkland (non-GE certified canola) and 11canola cultivars (L120, InVigor® 5440, PV 533 G, V22-1, L159, L252, PV 530 G, 74-44 BL, L150, L156H and 1022 RR) were received from Oilseeds Program of the Grain Research Laboratory of the Canadian Grain Commission. The 11 canola cultivars were used to assess consistency of the four endogenous reference sequences in terms of copy numbers measured. Breeder seeds of GT73 Conquest canola cultivar were received from Agricore Canada in 2001 and kept in a fridge. Certified reference material of non-GE Eagle canola cultivar was received from AOCS (Urbana, Illinois).

DNA extraction

The four DNA extraction kits used were: Fast ID Genomic DNA Extraction Kit (Genetic ID, Inc., Fairfield, IA), GM Quicker 2 (Nippon Gene Co., Ltd, Tokyo, Japan), NucleoSpin Food (Macherey-Nagel, D-Mark Biosciences, Toronto, ON, Canada), and DNeasy mericon Food (Qiagen Sciences, LLC, Louisville, KY). Manufacturer’s instructions were followed to extract DNA from 0.2 g samples. The amount of DNA was determined with PicoGreen Assay kit (Molecular Probes, Eugene, Oregon). Quantification assays were performed in 96-well fluorescence microtiter plates (Thermo Fisher Scientific, Waltham, MA, USA) and fluorescence was measured on Spectramax M5 Plate Reader (MDS Analytical Technologies, Toronto, Canada). Non-GE and GE canola DNA samples were mixed to prepare 0.01, 0.1 and 1% GE samples.

Digital PCR

QX200 ddPCR system (Bio-Rad, Pleasanton, CA) was used for the experiments. The primer and probe DNA sequences for the target and reference genes and the concentrations used for ddPCR and real-time quantitative PCR are listed in Table 1. Duplex ddPCR (mixing of target and reference primers and probes in the same reaction) was used to generate target and reference droplets at the same time. 12.5 μL Bio-Rad ddPCR Supermix for Probes (no dUTP) and either 10 ng or 100 ng DNA was used for the ddPCR assays. DG-32 cartridge for automated droplet generator (Cat. No. 186-4108) was used to generate droplets in 25 μL volume. The generated droplets were transferred to semi-skirted 96 well Eppendorf plates (Cat. No. 12001925) and sealed with pierceable foil heat seal (Cat. No. 1814040). MJ Thermal Cycler (PTC 200) was used for PCR amplification of the generated droplets. The thermal cycling conditions used were: 95 °C for 10 min initially, and then 50 cycles at 95 °C for 15 s and 60 °C for 1 min and finally 98 °C for 10 min. A ramp rate of 0.6 °C/sec was used between the cycling steps and a ramp rate of 0.3 °C/s was used at the last step to keep the reaction at 15 °C. The droplets were counted using the droplet reader of QX200 system. QuantaSoft version 1.7.4.0917 and automatic threshold were used for ddPCR data analysis. Bonferroni multiple-mean comparison method was used to compare selected means (https://www.spcforexcel.com/knowledge/comparing-processes/bonferronis-method).

Table 1

DNA sequences of primers and probes and concentrations used for ddPCR and real-time PCR.

Name of target event or reference gene	Sequences (5’ to 3’)	concentration used for ddPCR (μM)	Concentration used for real-time PCR (μM)
GT73	F-CCA TAT TGA CCA TCA TAC TCA TTG CT	0.4	0.15
	R-GCT TAT ACG AAG GCA AGA AAA GGA	0.4	0.15
	P-FAM-TTC CCG GAC ATG AAG ATC ATC CTC CTT-BHQ1	0.2
	P-FAM-TTC CCG GAC ATG AAG ATC ATC CTC CTT-TAMRA		0.05
HCN92	F-GTT GCG GTT CTG TCA GTT CC	0.4
	R-CGA CCG GCG CTG ATA TAT GA	0.4
	P-FAM-TCC CGC GTC ATC GGC GG-BHQ1	0.2
OXY235	F-GAT AGA TGG TGG TGT GAG TCT TGT	0.4	0.3
	R-CCT AAC TTT TGG TGT GAT GAT GCT	0.4	0.3
	P-FAM-TGC CAT CAG CTG ACA CGC CGT GC-BHQ1	0.2
	P-FAM-TGC CAT CAG CTG ACA CGC CGT GC-TAMRA		0.15
CruA	F-GGC CAG GGT TTC CGT GAT	0.2	0.2
	R-CCG TCG TTG TAG AAC CAT TGG	0.2	0.2
	P-HEX-AGT CCT TAT GTG CTC CAC TTT CTG GTG CA-BHQ1	0.2
	P-VIC-AGT CCT TAT GTG CTC CAC TTT CTG GTG CA-TAMRA		0.2
FatA(A)	F-ACA GAT GAA GTT CGG GAC GAG TAC	0.3	0.3
	R-CAG GTT GAG ATC CAC ATG CTT AAA TAT	0.9	0.9
	P-HEX- AAG AAG AAT CAT CAT GCT TC-BHQ1	0.15	0.15
HMG	F-GGT CGT CCT CCT AAG GCG AAA G	0.2	0.2
	R-CTT CTT CGG CGG TCG TCC AC	0.2	0.2
	P-VIC-CGG AGC CAC TCG GTG CCG CAA CTT-BHQ1	0.2	0.2
Ccf	F-ATT GGG CTA CAC CGG GAT GTG T	0.2
	R-GCT TCC GTG ATA TGC ACC AGA AAG	0.2
	P-HEX-CGA TGG TGT CCC CAG TCC TTA TGT GCT C-BHQ1	0.2

Real-time PCR

Real-time PCR was carried out using ABI 7500 PCR instrument (Thermo Fisher Scientific). 100 ng DNA (5 μL of 20 ng/ μL DNA solution), primer/probe concentrations provided in Table 1 and 1X (12.5 μL) TaqMan Universal Master Mix II with UNG (Applied Biosystems) were used for the assay. The thermal profile used for the real-time PCR was initial hold for 2 min at 50 °C, 10 min at 95 °C, and then 40 cycles of 15 s at 95 °C and 1 min annealing at 60 °C. The variation in thermal profile for real-time PCR and ddPCR was because of validation carried out at different times.

Results and discussion

Copy number variation among the commonly used canola reference genes

The copy numbers of HMG-I/Y and FatA(A) endogenous reference genes were significantly lower than the copy numbers of CruA and Ccf (Table 2). There was variability among the three DNA sources for each of the four reference targets. The copy numbers measured for Parkland cultivar were significantly higher than that of Legend and Eagle cultivars for HMG-I/Y and FatA(A) (Table 2). The haploid genome size of B. napus (canola) is estimated to range between 1129 and 1235 [10]. Based on 1.129 and 1.235 pg weight per haploid genome, and if 100 ng canola genomic DNA is used for PCR, there will be approximately 80,000–88,000 haploid genomic DNA copies. The use of 100 ng genomic DNA for PCR has also been reported to correspond to approximately 87,000 haploid copies of the B. napus genome [9]. For Legend and Eagle non-GM canola DNA, the number of copies for HMG-I/Y and FatA(A) varied from 6000 to 7000 for 10 ng DNA (Table 2). On the other hand, the number of copies for CruA and Ccf varied from 11,800 to 13,500 for 10 ng DNA. HMG-I/Y and CruA have been reported to be single and two copy genes, respectively [5]. A drawback reported for HMG-I/Y was that it may not be stable to be used as a reference gene for canola [4].

Table 2

Copy number variation of four endogenous references among three non-GE canola cultivars.

Canola cultivar (DNA source)	HMG-I/Y copy numbersX,Y	FatA(A) Copy numbersX,Y	CruA copy numbersX,Y	Ccf copy numbersX,Y	Comparison of reference means for each cultivarZ
Legend	6350 ± 217c	6083 ± 295b	11885 ± 462b	11852 ± 186b	[b, b, a, a]
Eagle	7283 ± 374b	6145 ± 192b	13369 ± 356a	13584 ± 460a	[b, b, a, a]
Parkland	8444 ± 149a	11764 ± 284a	13084 ± 552a	13397 ± 447a	[c, b, a, a]

10 ng DNA extracted with Fast ID DNA extraction method was used for ddPCR.

X Average of four ddPCR measurements plus minus standard deviation.

Y For each reference gene and copy numbers for the three cultivars, means assigned the same letter vertically are not significantly different (α = 0.05).

Z For each cultivar, reference means assigned the same letter horizontally are not significantly different (α = 0.05). For example, for Legend, means of HMG-I/Y, FatA(A), CruA and Ccf have the letters ‘b’, ‘b’, ‘a’ and ‘a’, respectively – means of CruA and Ccf were significantly higher than that of HMG-I/Y and FatA(A).

Effect of DNA extraction methods and DNA quality on copy number variation of reference genes

DNA samples extracted with four kits were used for ddPCR to determine the effect on the assessment of reference gene copy numbers. The reference gene copy numbers for Fast ID extracted DNA were significantly higher than DNA extracted with the other three extraction kits for all four reference genes (Table 3). A similar trend was observed for both Legend and Eagle DNA. There was relatively less variability of copy numbers for DNA extracted with GM Quicker II, DNeasy® mericon Food and NucleoSpin Food. HMG-I/Y and FatA(A) had lower copy numbers than CruA and Ccf (Table 3). The results indicated that different DNA extraction methods may cause variation in copy numbers, and thus it is important not to use DNA extracted with different methods for the same experiment.

Table 3

Effect of four DNA extraction methods on the assessment of reference gene copy numbers.

DNA extraction method and cultivar	HMG-I/Y copy numbers	FatA(A) copy numbers	CruA copy numbers	Ccf copy numbers
FID – Legend	6841 ± 149a	6463 ± 43a	13741 ± 97a	13704 ± 304a
GMQ2 – Legend	4746 ± 234b	4662 ± 60bc	9641 ± 288b	9898 ± 252b
DMF – Legend	4673 ± 116b	4801 ± 31b	9221 ± 106b	9221 ± 276c
NSF – Legend	4738 ± 54b	4548 ± 127c	9198 ± 289b	9684 ± 137bc
FID – Eagle	7018 ± 202a	6550 ± 137a	13943 ± 578a	14417 ± 249a
GMQ2 – Eagle	4866 ± 89b	4610 ± 78b	9730 ± 193bc	9989 ± 246b
DMF – Eagle	4487 ± 228b	4526 ± 5b	8649 ± 296c	9009 ± 277c
NSF – Eagle	5018 ± 525b	4641 ± 93b	9963 ± 222bc	10560 ± 392b

FID = Fast ID DNA extraction kit; GMQ2 = GM Quicker II DNA extraction kit; DMF = DNeasy® mericon Food kit; NSF = NucleoSpin Food kit. Average of three ddPCR measurements plus minus standard deviation. For each reference gene/cultivar and four DNA extraction methods (vertically), means assigned the same letter are not significantly different (α = 0.05). 10 ng DNA was used for ddPCR.

The copy number of three reference genes was determined for non-GE Legend canola DNA extracted at various times (Table 4). The Legend DNA extracted in 2015 was more sheared than the three DNA extracted in 2017 (Fig. 1). We have observed that Fast ID extraction method has a tendency to shear DNA of canola samples. The highly sheared Legend DNA had at least twice the number of copies for HMG-I/Y, FatA(A) and CruA reference genes and was significantly higher than the copy number of the three Legend DNA that were relatively less sheared. Thus, absolute quantification of GE DNA with digital PCR can be affected if the DNA is sheared. We do not have a concrete explanation for this. It is possible that some of the sheared DNA is denatured. According to Huggett et al. 2013 [11], double stranded DNA (dsDNA) molecule occupies one partition. If dsDNA is denatured, two single strands are formed that occupy two partitions. The Bio-Rad manual recommends shearing of DNA using restriction enzymes. Shearing of DNA reduces sample viscosity and improves template accessibility that enhances assay performance [12].

Table 4

Copy number variation for different Legend DNA sources.

DNA source	HMG-I/Y copy numbers	FatA(A) copy numbers	CruA copy numbers
Legend-1	6670 ± 137b	6614 ± 217b	12448 ± 570b
Legend-2	6547 ± 175b	6104 ± 177bc	12580 ± 172b
Legend-3	6271 ± 140b	5989 ± 75c	11698 ± 123b
Legend-4	16212 ± 366a	13846 ± 355a	30819 ± 706a

DNA of Legend 1–4 were extracted in April 2017, March 2017, January 2017 and October 2015, respectively. Fast ID DNA extraction method was used. Fig. 1 shows the quality of DNA extracted from the four Legend sources. Average number of copies for four ddPCR measurements plus minus standard deviation. For each reference gene and the four DNA sources, means assigned the same letter vertically are not significantly different (α = 0.05). 10 ng DNA was used for ddPCR.

Fig. 1

Agarose gel-electrophoresis (1.2%) of Legend canola DNA extracted at different times. Fast ID Genomic DNA Extraction kit was used. M = Low DNA Mass Ladder (2 μL of 117.5 ng/μL was loaded on the gel); 1–4 indicate 400 ng DNA extracted at four different times (April, March, January 2017 and October 2015, respectively).

Stability of reference copy numbers among canola cultivars

DNA was extracted from eleven canola cultivars in order to assess copy number variation among the four endogenous references used (Table 5). Most variability among the canola cultivars was observed for HMG-I/Y reference gene (Fig. 2). Our observation agrees with that of Wu et al. [4] in that HMG-I/Y may not produce consistent results among different canola cultivars. Instability of HMG-I/Y was attributed to single nucleotide polymorphism (SNP) sites in the primer DNA sequence [4]. On the other hand, consistent and repeatable results were obtained for FatA(A) reference gene. Thus, FatA(A) primer/probe DNA sequences seem to be suitable for absolute quantification of GE canola events by ddPCR. Our observation confirms the work of Henderson et al. [9] who reported FatA(A) to be a suitable reference gene for real-time quantitative PCR analysis of GE canola events. An example of ddPCR amplitude plot for CruA and FatA(A) endogenous reference genes is provided in Fig. 3. There was a clear separation of positive and negative clusters.

Table 5

Copy number variation of reference targets among 11 canola cultivars.

Cultivar name	HMG-I/Y copy numbers	FatA(A) copy numbers	CruA copy numbers	Ccf copy numbers
L1120	6823 ± 81	5746 ± 213	11955 ± 288	12507 ± 373
InVigorR 5440	3577 ± 113	5895 ± 227	12138 ± 355	12511 ± 420
PV 533 G	5515 ± 216	5051 ± 141	9842 ± 310	10182 ± 282
V22-1	3131 ± 25	5382 ± 295	11110 ± 116	11469 ± 589
L159	3484 ± 57	5296 ± 308	10613 ± 133	10786 ± 330
L252	3277 ± 191	5580 ± 36	11973 ± 559	12202 ± 155
PV 530 G	3099 ± 142	5229 ± 149	10540 ± 151	10789 ± 295
74-44 BL	6037 ± 154	5420 ± 171	10904 ± 267	10845 ± 379
L150	2951 ± 41	5292 ± 42	10523 ± 137	10714 ± 371
L156H	5845 ± 153	5104 ± 194	10431 ± 282	10697 ± 569
1022 RR	3407 ± 237	5169 ± 120	10826 ± 268	10848 ± 141
Overall mean	4286 ± 1405	5379 ± 303	10987 ± 754	11232 ± 848

Average number of copies for three ddPCR measurements plus minus standard deviation. 10 ng DNA extracted with Fast ID DNA extraction kit was used for ddPCR.

Fig. 2

Average number of copies obtained for reference genes used (average of 11 canola cultivars, n = 33) plus minus standard deviation. 10 ng DNA was used for ddPCR.

Fig. 3

Examples of amplitude plots for ddPCR assays. The amplitude plots are for Legend canola cultivar and FatA(A) (single copy) and CruA (two copy) reference genes. Green = positive droplets; black = negative droplets.

Droplet digital PCR and real-time quantitative PCR results for GE canola samples

Accuracy of ddPCR for absolute quantification of 0.01%, 0.1% and 1% OXY235, GT73 and HCN92 canola DNA samples was assessed using three reference genes and two DNA extraction methods (Table 6). Non-GE template controls were included for analysis of spiked GE samples for each event, and no signals were observed for the negative template controls. Overall, the percentage values obtained for CruA were much lower than that of HMG-I/Y and FatA(A). Dividing OXY235, GT73 and HCN92 target copy numbers by CruA reference copy numbers resulted in lower than expected percentage values as more copies were generated for the two copy CruA reference gene. The results indicated that if a reference gene with two copies is used for ddPCR assessment of GE canola events, erroneous quantification results may be achieved unless a correction is made. There were also some variations between the two DNA extraction kits (Table 6). Some of the results for 0.01% GE DNA samples were inconsistent as the number of target copies were either few or no target copies were obtained for some replications. For a comparison, we carried out real-time quantitative PCR experiments to assess the percentage values obtained for OXY235 and GT73 GE canola samples using CruA, FatA(A) and HMG-I/Y reference genes. Overall, the expected RT-PCR results were achieved for both CruA (two copy) and FatA(A)/HMG-I/Y (single copy) reference genes (Table 7). The observation indicated that variation in copy number for reference genes did not affect real-time PCR results. On the other hand, ddPCR results were affected as discussed above.

Table 6

Droplet digital PCR results obtained for GE canola samples using two DNA extraction methods and three reference genes.

Event	DNA extraction method	Reference gene	0.01% GE sample	0.1% GE sample	1% GE sample
OXY235	FID	HMG-I/Y	0.005 ± 0.003	0.08 ± 0.02	0.68 ± 0.03
		FatA(A)	0.008 ± 0.002	0.066 ± 0.006	0.68 ± 0.06
		CruA	0.004 ± 0.002	0.029 ± 0.01	0.349 ± 0.04
	DMF	HMG-I/Y	0.006 ± 0.002	0.104 ± 0.03	1.00 ± 0.09
		FatA(A)	0.009 ± 0.008	0.10 ± 0.01	0.90 ± 0.02
		CruA	0.006 ± 0.003	0.05 ± 0.019	0.419 ± 0.01
GT73	FID	HMG-I/Y	0.007 ± 0.004	0.096 ± 0.02	0.866 ± 0.07
		FatA(A)	0.01 ± 0.004	0.084 ± 0.012	0.931 ± 0.05
		CruA	0.005 ± 0.003	0.043 ± 0.007	0.446 ± 0.02
	DMF	HMG-I/Y	0.018 ± 0.007	0.086 ± 0.01	0.838 ± 0.03
		FatA(A)	0.016 ± 0.002	0.099 ± 0.01	0.827 ± 0.05
		CruA	0.007 ± 0.001	0.048 ± 0.03	0.429 ± 0.03
HCN92	FID	HMG-I/Y	0.009 ± 0.001	0.115 ± 0.009	1.08 ± 0.05
		FatA(A)	0.011 ± 0.003	0.100 ± 0.007	1.10 ± 0.1
		CruA	0.07 ± 0.003	0.052 ± 0.002	0.520 ± 0.24
	DMF	HMG-I/Y	0.006 ± 0.001	0.064 ± 0.04	0.684 ± 0.07
		FatA(A)	0.005 ± 0.006	0.07 ± 0.015	0.728 ± 0.05
		CruA	0.005 ± 0.005	0.035 ± 0.008	0.348 ± 0.03

FID = Fast ID DNA extraction kit; DMF, DNeasy® mericon Food kit. Average of three ddPCR measurements plus minus standard deviation. 100 ng DNA was used for ddPCR.

Table 7

Example of real-time quantitative PCR results obtained for DNA extracted from GE canola samples using Mericon DNA extraction kit.

Canola event	Reference gene	0.01% GE sample	0.1% GE sample	1% GE sample
OXY235	HMG-I/Y	0.011 ± 0.01	0.11 ± 0.015	0.82 ± 0.06
	FatA(A)	0.012 ± 0.00	0.08 ± 0.01	0.75 ± 0.01
	CruA	0.06 ± 0.00	0.08 ± 0.01	0.75 ± 0.06
GT73	HMG-I/Y	0.009 ± 0.00	0.08 ± 0.02	0.66 ± 0.01
	FatA(A)	0.008 ± 0.00	0.07 ± 0.01	0.58 ± 0.01
	CruA	0.013 ± 0.00	0.084 ± 0.01	0.77 ± 0.02

Average of three qPCR measurements plus minus standard deviation. 100 ng DNA extracted with Fast ID DNA extraction kit was used for real-time quantitative PCR.

Conclusions

Our results indicated the importance of using appropriate reference gene for absolute quantification of GE canola events by ddPCR. HMG-I/Y and FatA(A) can be used as a reference gene for absolute quantification of GE canola events by ddPCR. However, FatA(A) is more consistent than HMG-I/Y and will be the preferred reference gene for ddPCR assessment of GE canola events. The expected GE percentage values will be underestimated if two copy reference gene (CruA and Ccf) is used for ddPCR. Real-time quantitative PCR analysis of GE samples was not affected by the use of single vs. two copy reference genes. It is also important to pay attention to the DNA extraction method and quality of DNA used for ddPCR.

Conflict of interest

The authors declare that there are no conflicts of interest