All animal procedures were carried out after approval of the Institutional Animal Care and Use Committee (IACUC) of Optipharm, Inc., Life Science Institute (IACUC approval No. OPT-140103-1).

TALE sequences and TALEN target site was identified by using CHOPCHOP, a web tool for selecting target sites for TALEN-directed mutagenesis (http://chopchop.cbu.uib.no) [12]. Exon 9 of porcine GGTA1 gene (GenBank accession no. AH010595.2) was chosen for TALEN construct with 5′-GAG GAG TTC TTA A-3′ for left TALE, 5′-GGC CAC AAA GTC ATC TTT TAC A-3′ for right TALE. Constructed TALEN was transfected into primary ear skin fibroblasts of Yucatan male miniature pig at 3-month old using nucleofector, Amaxa-4D (Lonza, Köln, Cologne, Germany). After transfection, 1×10⁴ transfected cells were seeded in a 100-mm tissue culture plate (Corning, NY, USA) for single colony culture. Genomic DNA was extracted from each colony using DNeasy Blood and Tissue kit (Qiagen, Redwood City, CA, USA). Target sequence of TALEN in Exon 9 of GGTA1 was amplified with 5′-AGA ATC ACC AGT CAG GTA AGC CAC TCC-3′ (forward primer) and 5′-TTG GAA GAC CTG ATC CAC GTC CAT GCA G-3′ (reverse primer), and the amplicon was subject to Sanger sequencing to analyze the mutagenesis in target site.

Oocyte preparation, nuclear transfer, embryo culture and embryo transfer were performed according to Choi et al [13]. Briefly, Cumulus-oocyte complexes were collected from prepubertal gilts at a local abattoir and maturated in TCM 199 medium (Gibco, Waltham,, MA, USA) supplemented with: 0.5 μg/mL luteinizing hormone, 0.5 μg/mL follicle stimulating hormone, 10 ng/mL epidermal growth factor, plus 0.1% polyvinyl alcohol, 3.05 mM D-glucose, 0.91 mM sodium pyruvate, 0.57 mM cysteine, 75 μg/mL penicillin G, and 50 μg/mL streptomycin) at 39°C with 5% CO₂ air condition. After maturation, oocytes were freed from cumulus and enucleated by aspiration. TALEN-treated donor cells were cultured for 3 days in Dulbecco’s modified eagle medium supplemented with 0.5% fetal bovine serum, and a single donor cell was placed in the perivitelline space of an oocyte. Fusion and activation were induced simultaneously with two successive DC pulses of 1.1 kV/cm for 30 μs (BTX, Holliston, MA, USA). Reconstructed embryos were allowed for 1 or 2 days of additional culture in North Carolina State University (NCSU)-23 medium supplemented with 0.4% bovine serum albumin, and surgically transferred to the one lateral oviduct of a naturally cycling sow on the first day of standing estrus. Pregnancy status was monitored using an ultrasound scanner (Medison Co., Pangyo, Korea). F2 piglets of GGTA1 knockout pigs were produced by backcross between TALEN-mediated founder and F1 piglets. Ear skin fibroblasts were harvested from F2 piglets and wild type strain controls. The cells were stained with 1 μL of Alexa 488-conjugated isolectin GS-IB4 against (Invitrogen, Waltham, MA, USA) in 100 μL of Dulbecco’s phosphate-buffered saline (Gibco, USA) for 1 h at 4°C to detect alpha-gal epitope. After staining and washing, the cells were subject to flow cytometry analysis using FACSCallibur (Becton, Dickinson and Company, Franklin Lakes, NJ, USA).

Genomic DNA was extracted from ear skin of GGTA1 knockout pigs (CB1, CB3, CJ1) and wild type strain controls (WT1-3) using DNeasy Blood & Tissue Kit (Qiagen, USA). One microgram of genomic DNA sample was subjected to library preparation. Whole genome sequencing libraries with short inserts of 350 to 450 bp for paired-end reads were prepared using Truseq DNA sample prep kit (Illumina, San Diego, CA, USA) according to the manufacturer’s instruction. Whole genomic DNA was sequenced using Illumina HiSeq 2500 with an adaptation of the pairwise end-sequencing strategy. Post-sequencing analysis was performed using Theragen Etex’ workflow (www.theragenetex.com). Briefly, sequencing reads were aligned with Burrows-Wheeler aligner [14] with default parameters to The Swine Genome Sequencing Consortium Sscrofa 10.2 (https://doi.org/www.ncbi.nlm.nih.gov/assembly/GCF_000003025.5/). Alignments were sorted, merged and de-duplicated with SAMtools [15]. Genome Analysis Tool Kit was applied for variant calling with default parameters [16]. Variants were annotated by SnpEff [17] and depth of each base was calculated by BEDTools [18]. Genome-wide distribution of variant calls was visualized by using Circa (http://omgenomics.com/circa/).

Total RNA was extracted from whole blood of GGTA1 knockout pigs (CB1, CB3, CJ1) and wild type strain control (NC) using QIAamp RNA Blood Mini Kit (Qiagen, USA). RNA concentration and purity were determined using the Agilent 2100 BioAnalyzer (Agilent, Palo Alto, CA, USA), and RNA samples with RNA integrity number ≥7.0 were subject to library preparation. RNA sequencing libraries were prepared using 2 μg of total RNA as paired-end reads with a length of 100 bases using the TruSeq RNA Sample Preparation Kit (Illumina, USA). After quantification and qualification of the libraries using the Agilent 2100 BioAnalyzer (Agilent, USA) and the KAPA library quantification kit (Kapa Biosystems, Wilmington, MA, USA), flow cells were sequenced as paired-end reads (2×100 bp) using the Illumina HiSeq 2500 (Illumina, USA). Read filtering process and differentially expressed gene analysis were performed using Theragen Etex’ workflow (www.theragenetex.com). Briefly, after filtering out the low-quality reads, Trinity method was applied to assemble short reads [19], CD-HIT-EST to select the final unigenes from the assembled transcripts with identities of 94% [20]. Incomplete and redundant assembled transcripts were removed by selecting representative sequences [21]. Gene expression level was calculated based on fragments per kilobase of exon per million mapped reads (FPKM) using Cufflinks from The Swine Genome Sequencing Consortium Sscrofa 10.2. Differentially expressed genes were identified at a fold change cutoff of 2. Multiplex literature mining was performed using PubMatrix [22] with pair-wise comparisons of differentially expressed genes and modifier terms including nuclear transfer, lung, alveologenesis, and pig. False-positive findings were removed by manual curation.

Heterozygous GGTA1 knockout pig ear skin fibroblasts were generated by homologous recombination at exon 4 of the gene, and used to generate GGTA1 knockout pigs by nuclear transfer (piglet ID = CJ1) [13]. The same pig ear skin fibroblast clones were subjected to TALEN treatment for mutation of exon 9 (Figure 1A). After verification of genetic interruption and reduced expression of GGTA1 protein in the TALEN-mediated pig ear skin fibroblasts, the cells were used as donor cells for nuclear transfer, and reconstructed embryos were transferred to recipient miniature pigs. After full-term development, TALEN-mediated GGTA1 knockout piglets were successfully delivered (piglet IDs = CB1 and CB3). Sequencing analysis of genomic DNA isolated from the ear skin of the TALEN-mediated GGTA1 knockout piglets demonstrated that CB1 had various mutant genotypes, such as a single base (A) insertion at chr1:293631881 and a single base deletion (T) at chr1:293631880 in exon 9 of GGTA1 gene. In CB3, there was a 6-base deletion (ATA CTT) and a 47-base insertion (GTG GAG ACT TGG AAA TCC CCG TGA GTC AAA CCG CTA TCC ACG CCC ATT GAT GTA CTG C) at chr1: 293631888 – chr1:293631883 in exon 9 of the GGTA1 gene (Figure 1B). The inserted sequence in CB3 was not identified by the National Center for Biotechnology Information Basic Local Alignment Search Tool (https://doi.org/blast.ncbi.nlm.nih.gov/Blast.cgi). Small indels occurred between the left and right TALE sequences. To confirm germline transmission of the genetic traits of the TALEN-mediated GGTA1 knockout pigs, one F2 offspring of CB1 (piglet ID = CB1 F2) and three F2 offspring of CB3 (piglet IDs = CB F2-1 to 3) were produced by backcross. Flow cytometry analysis showed that the expression levels of GGTA1 protein in fibroblasts derived from the ear skin of F2 piglets were significantly lower than those of the wild-type strain control (Figure 1C). These results demonstrated that TALEN is an efficient and precise genome-editing method to modulate target gene activity, and that TALEN-mediated gene mutations are heritable to the next generation.

Although GGTA1 was successfully disrupted by TALEN, we could not exclude the possibility of off-target mutations in genetic elements not targeted by TALEN. To assess this, we performed whole-genome sequencing to determine the presence of indels and single nucleotide polymorphisms (SNPs) in the TALEN-mediated GGTA1 knockout pigs. We first compared the whole-genome sequence of donor cells to the pig reference genome Sscrofa 10.2, and found 802,662 indels and 3,878,732 SNPs throughout the chromosomes (Figure 2A). The unexpectedly high frequency of genomic mutations in the donor cells implied that the pig reference genome would be inappropriate for identification of TALEN-induced, true off-target mutations by sequence comparison. Therefore, we determined the whole-genome sequence of the donor cells for use as an internal control genome, to obtain a precise estimation of genomic mutations induced by TALEN in GGTA1 knockout pigs.

Mutation profiling using the internal control genome revealed that TALEN-mediated GGTA1 knockout pigs had few indels and SNPs, with numbers comparable to the homologous recombination-treated knockout pig and the wild-type strain controls (Figure 2B). The average sequencing depth was 230× for TALEN-treated pigs, 200× for homologous recombination-treated pigs, and 204× for wild-type controls, which confirmed the accuracy of the whole-genome sequencing analyses. The variant calling pipeline revealed that the TALEN-treated pig CB1 harbored one base insertion (T) at chr8: 33762665, one base deletion (T) at chr6:81281307, and eight SNPs, while CB3 had 70 indels and 192 SNPs. No off-target mutations induced by TALEN were predicted by CHOPCHOP [12] in the basis of Sscrofa 10.2. Compared to the TALEN-treated pigs, the homologous recombination-treated pig CJ1 had a single deletion (T) at chr8: 33762665 and 11 SNPs, while the wild-type controls had an average of 80 indels and 223 SNPs. The mutation rates detected in the TALEN-mediated GGTA1 knockout pigs were similar or higher than those in the homologous recombination-treated pig, but lower than the frequency of spontaneous germline variants observed in wild-type strain controls.

The genomic mutations were distributed evenly throughout the chromosomes in TALEN- and homologous recombination-treated pigs, and were comparable to those of germline variants in wild-type strain controls (Figure 2C). Among the 10 genomic mutations found in CB1, one indel and two SNPs were in intron regions of protein-coding genes. Of the 20 indels in CB3, two were found in the upstream region, one in the 3′ untranslated region (UTR), and 17 in the intron region. Out of 59 SNPs, six were located in the upstream region, one in the 3′ UTR, and 52 in the intron regions of protein-coding genes. Similarly, SNPs detected in CJ1 were located in the intron regions of protein-coding genes. None of the mutations in the TALEN- or homologous recombination-treated pigs were expected to have deleterious effects (Table 1, Supplementary Tables S1–6). Taken together, the genomic analysis of GGTA1 knockout pigs demonstrated that an internal control genome is important for precise detection of true genomic mutations due to genome editing, and that TALEN might not induce excessive off-target mutations, with frequencies comparable to those associated with homologous recombination.

As described above, TALEN-induced mutations in GGTA1 knockout pigs were not excessive compared to the internal control genome, and most of the mutations were located in non-coding regions of various genes; however, this does not exclude the possibility of unintended effects of TALEN-induced mutations on the physiological and genomic integrity of GGTA1 knockout pigs. To address this issue, we performed transcriptomic analysis of GGTA1 knockout pigs, and identified differentially expressed genes in those pigs by comparison with the wild-type strain control.

Expression values of an average of 56,038 transcripts from 22,852 genes were obtained from GGTA1 knockout and wild-type pigs (Supplementary Table S7). Correlation analysis of the expressed transcripts showed that the global gene expression profiles of individual knockout pigs and the wild-type pig were well-correlated (correlation coefficient of 0.744 to 0.857; Figure 3A). In addition, few genes were differentially expressed in GGTA1 knockout pigs compared to the wild-type pig (fold change ≥2). CB1 had 23 up- and seven down-regulated genes compared to the wild-type pig. Likewise, CB3 had 39 up- and six down-regulated genes. The numbers of differentially expressed genes in the TALEN-mediated GGTA1 knockout pigs were comparable to those in the homologous recombination-treated CJ1 pig, which had 35 up- and five down-regulated genes.

We compared the differentially expressed genes of the individual knockout pigs, and determined that three, six, and four genes were exclusively expressed in CB1, CB3, and CJ1, respectively. In contrast, 35 genes were commonly expressed in GGTA1 knockout pigs (Figure 3B), including Shroom family member 2 (SHROOM2), cysteine and glycine rich protein 3, and MARVEL domain-containing 2 (MARVELD2) (Figure 3B). When we compared the fold-changes of the genes commonly differentially expressed in GGTA1 knockout pigs, most showed similar directions of change (up- or down-regulation) compared to the wild-type pig (Figure 3C). Literature mining revealed that several differentially expressed genes in GGTA1 knockout pigs are associated with DNA methylation in gene regulation, including family with sequence similarity 50, member B (FAM50B) [23,24], cytochrome C, somatic [25], SHROOM2 [26], protein tyrosine phosphatase, non-receptor type 13 (PTPN13) [27], and LDL receptor-related protein 12 (LRP12) [28]. Other differentially expressed genes are known to be involved in respiratory diseases, such as lung cancer, including FAM50B [24], ATPase phospholipid transporting 11C [29], WWC family member 3 [30–32], PTPN13 [27], RAS protein activator-like 2 [33], tyrosine 3-monooxygenase/tryptophan 5-monooxygenase activation protein zeta [34,35], and LRP12 [28]. Dysregulated DNA methylation and lung pathophysiology are well known abnormalities in pigs cloned by somatic cell nuclear transfer [36]. Therefore, it appears that the differential gene expression in GGTA1 knockout pigs is caused by the nuclear transfer procedure itself, rather than off-target effects of TALEN.

To obtain insight into the association between TALEN-induced genomic mutations and gene expression, we generated multi-layer Circos plots with differentially expressed transcripts, SNPs, and indels. GGTA1 knockout pigs produced by TALEN had 1,873 differentially expressed transcripts covering seven annotated genes (CB1), and 3,722 transcripts with 16 annotated genes (CB3) (log2|fold change|≥2). In addition, CB1 had eight SNPs and two indels, while CB3 had 192 SNPs and 70 indels. All of the differentially expressed transcripts, SNPs, and indels from TALEN-mediated GGTA1 knockout pigs were distributed evenly throughout the genome (Figure 4A, B). Likewise, the GGTA1 knockout pig produced by homologous recombination had 3,142 differentially expressed transcripts covering 14 annotated genes, 11 SNPs, and one indel, which were distributed evenly throughout the genome (Figure 4C). Among the genomic mutations, one SNP (c.376 +331C>T) at chr3:33134141 in the intron region of ENSSSCT 00000008657 was common to both CB3 and CJ1, and one insertion (c.820-502_820-501insT) at chr6:81281307 in the intron region of syndecan-3 (SDC3) was common to CB1 and CB3. There were no other common genomic mutations among the GGTA1 knockout pigs.

Next, we investigated the possible effects of genomic mutations on gene expression by comparing the expression levels of protein-coding genes harboring SNPs or indels with those of intact genes (Figure 4D). As described above, CB1 and CB3 shared a single base insertion in the intron region of SDC3; their expression levels were comparable to that of SDC3 in CJ1, which had no mutation in that region. Likewise, CB3 had indels in non-coding regions of Staufen double-stranded RNA-binding protein 2, NEDD4-binding protein 2, signal-induced proliferation-associated 1-like 2, zinc finger HIT-type containing 3, collagen type XIX alpha 1 chain, SUMO-specific peptidase 6, and transcription factor 7-like 2; the expression levels of these genes were comparable to those in CB1 and CJ1, which had no genomic mutations in these genetic elements. Furthermore, no differences in expression levels were observed among genes harboring SNPs in non-coding regions in CB3 compared to CB1 and CJ1, including MARVELD2, X-ray repair cross-complementing 2, and ribosomal protein L31. These results imply that TALEN-induced off-target mutations in GGTA1 knockout pigs had no discernable effects on RNA transcript abundance.