Transcriptome Comparison of Human Neurons Generated Using Induced Pluripotent Stem Cells Derived from Dental Pulp and Skin Fibroblasts

Induced pluripotent stem cell (iPSC) technology is providing an opportunity to study neuropsychiatric disorders through the capacity to grow patient-specific neurons in vitro. Skin fibroblasts obtained by biopsy have been the most reliable source of cells for reprogramming. However, using other somatic cells obtained by less invasive means would be ideal, especially in children with autism spectrum disorders (ASD) and other neurodevelopmental conditions. In addition to fibroblasts, iPSCs have been developed from cord blood, lymphocytes, hair keratinocytes, and dental pulp from deciduous teeth. Of these, dental pulp would be a good source for neurodevelopmental disorders in children because obtaining material is non-invasive. We investigated its suitability for disease modeling by carrying out gene expression profiling, using RNA-seq, on differentiated neurons derived from iPSCs made from dental pulp extracted from deciduous teeth (T-iPSCs) and fibroblasts (F-iPSCs). This is the first RNA-seq analysis comparing gene expression profiles in neurons derived from iPSCs made from different somatic cells. For the most part, gene expression profiles were quite similar with only 329 genes showing differential expression at a nominally significant p-value (p<0.05), of which 63 remained significant after correcting for genome-wide analysis (FDR <0.05). The most striking difference was the lower level of expression detected for numerous members of the all four HOX gene families in neurons derived from T-iPSCs. In addition, an increased level of expression was seen for several transcription factors expressed in the developing forebrain (FOXP2, OTX1, and LHX2, for example). Overall, pathway analysis revealed that differentially expressed genes that showed higher levels of expression in neurons derived from T-iPSCs were enriched for genes implicated in schizophrenia (SZ). The findings suggest that neurons derived from T-iPSCs are suitable for disease-modeling neuropsychiatric disorder and may have some advantages over those derived from F-iPSCs.


Introduction
We and other groups are using induced pluripotent stem cells (iPSCs) for in vitro disease modeling in a variety of neuropsychiatric disorders, including schizophrenia (SZ) and autism spectrum disorders (ASD) [1][2][3][4][5][6][7][8][9][10]. In addition to their utility for disease modeling in terms of identifying patient vs control differences in gene expression, morphology, synaptic architecture, and neuronal function, iPSCs can also be used to study human neurogenesis in vitro, which is particularly relevant to SZ and ASD considering that both have a neurodevelopmental basis. A variety of cell types have been used for iPSC reprogramming, but fibroblasts obtained from skin biopsy samples have been the mainstay for neuropsychiatric disorders so far. This presents a potential obstacle for modeling genetically-based childhood disorders. Although iPSCs have been developed from children with ASD and other developmental problems using fibroblasts and more recently, peripheral blood, it is somewhat problematic because such children often fear medical procedures, even routine phlebotomy [3][4][5]8]. Thus, alternative sources of cells for iPSC reprogramming obtained by non-invasive means would be useful. Also, because iPSCs may retain some cell-of-origin epigenetic marks, testing  differentiating neurons derived from iPSCs generated from various somatic cells to assess their utility for modeling neuropsychiatric disorders is important. One potential source of somatic cells is dental pulp derived from deciduous teeth. iPSCs derived from dental pulp (T-iPSCs) display typical molecular and cellular features of pluripotency, and have been shown to differentiate into neurons [11][12][13]. However, a more detailed molecular profile is needed to assess the similarity of iPSC-derived neurons from different somatic tissues/cells and their potential use for modeling neuropsychiatric disorders. Consequently, we have carried out an extensive gene expression profiling analysis of neurons derived from T-iPSCs using whole transcriptome profiling (RNA-Seq) and compared that with neurons derived iPSCs made from fibroblasts (F-iPSCs). Our gene expression profiling studies show a high degree of correlation for the two sources of neurons. However, there are subtle differences that might influence the decision to use T-iPSCs or F-iPSCs for some neuropsychiatric disorders.

Development of iPSCs from Dental Pulp and Skin Fibroblasts
The study was approved by the Albert Einstein College of Medicine Institutional Review Board. Written informed consent was obtained for subjects undergoing a skin biopsy, which was carried out by a board-certified dermatologist. For the tooth sample, signed written assent was provided by the subject, who was 12 years old at the time; the assent was countersigned by a parent. The tooth sample was lost naturally and not extracted. Consent for the skin biopsy samples was obtained by a senior member of the research team (the corresponding author). Assent for the tooth sample was obtained by a senior level associate (Ph.D. level) of one of the co-authors (JJF). T-iPSC lines (TIPS4 and TIPS4-C5) were generated from dental pulp cells harvested from a molar tooth that was naturally shed by the subject, a 12 year old healthy Caucasian male. The method for collecting deciduous teeth, extracting dental pulp, and growing these cells in culture is described in greater detail in Supplemental Methods (Text S1). Fibroblasts were obtained from skin biopsies performed in consenting adults by a board-certified dermatologist. The detailed procedure for growing fibroblasts in preparation for reprogramming into iPSCs is also in the Text S1. The F-iPSC lines referred to throughout the paper as F-iPSC1 and F-iPSC2 were derived from a 30 year old female and a 58 year old male, respectively. The development of these lines was previously described [1,2,14].
Immunocytochemistry was carried as previously described [19,20]. A list of the antibodies used in the study is shown in Text S1.

Neuronal Differentiation
Neurons were derived from neural progenitor cells (NPCs) as described by Marchetto et al. with slight modifications [1,9]. A detailed description of the protocol is in Text S1.

RNA-Seq
RNA-seq was carried out on iPSCs, NPCs and day 14 neurons derived from TIPS4 and TIPS4-C5, and day 14 neurons from F-iPSC1 and F-iPSC2. Total RNA was isolated from cells using the miRNeasy Kit (Qiagen) according to the manufacturer's protocol. An additional DNAse1 digestion step was performed to ensure that the samples were not contaminated with genomic DNA. RNA purity was assessed using the Agilant 2100 Bioanalyzer (Beijing Genomics Institute). Each RNA sample had an A260:A280 ratio above 1.8, a RIN.9, and an A260:A230 ratio above 2.2. Briefly, total RNA was converted to cDNA using oligo dT, which was then used for Illumina sequencing library preparation. Paired end RNA-seq was carried on an Illumina HiSeq 2000. We obtained 90-bp mate-paired reads from DNA fragments of with an average size of 250-bp (standard deviation for the distribution of inner distances between mate pairs is approximately 50 bp). RNA-Seq reads were aligned to the human genome (GRCh37/hg19) using the software TopHat (version 2.0.8) [21]. We counted the number of fragments mapped to each gene annotated in the GENCODE database (version 15) [22]. The category of transcripts is described at http://vega.sanger.ac.uk/info/about/ gene_and_transcript_types.html. Transcript abundances were  measured in FPKM (fragments per kilobase of exon per million fragments mapped). We used DESeq (an R package developed by Anders and Huber) to evaluate differential expression from count data [23]. We used DESeq (an R package developed by Anders and Huber) to evaluate differential expression from count data {{2762 Anders,S. 2010}}. Specifically, DESeq models the variance in fragment counts across replicates using the negative binomial distribution and tests whether, for a given gene, the change in expression strength between the two experimental conditions is significantly large as compared to the variation within each replicate group. In the end, only genes with average FPKMs larger than 1 across samples were considered for differential expression. The number of reads obtained from the RNA-seq runs for each sample and the fraction that could be aligned to the human genome was consistent across samples (Table S1). In addition, the correlation coefficients were very high for biological replicates (Table S1). Sequence data can be accessed at the NCBI's (National Center for Biotechnology Information) Gene Expression Omnibus (accession number GSE43143).

Reverse Transcribed PCR (RT-PCR) and Quantitative Realtime PCR (qPCR)
Reverse transcribed PCR (RT-PCR) was performed using a OneStep RT-PCR Kit (Qiagen, Valencia, CA) according to the manufacturer's instructions. The cDNA was generated using an iScriptTM cDNA synthesis Kit (Bio-RAD, Hercules, CA) and subsequently used as a template for quantitative qPCR, which was carried out with an ABI 7900HT Real-Time PCR System instrument (Applied Biosystems, Foster City, CA). Each reaction consisted of cDNA, primers, and SYBR Green PCR Master Mix (Applied Biosystems, Foster City, CA) in an 8 ml volume (primers used in this study are shown in Text S1). Relative changes in gene expression were calculated using the 2 2DDCt method with b2-microglobulin (b2M) as a reference gene. The primers used in this study, as well as technical details are in Text S1.

Electrophysiology
Whole cell recordings were made using a Multiclamp 700B (Molecular Devices, Sunnyvale, CA) in 35 day old differentiated cultures. Neuronal-like cells characterized by extended processes were chosen for recording. Low-resistance pipettes (3-5 MOhm) contained: 135 mM K gluconate, 6 mM NaCl, 10 mM HEPES, 1 mM EGTA, 0.5 mM CaCl2, 10 mM Glucose, 2 mM MgATP, 0.3 mM NaGTP pH 7.2. Cells were perfused with an extracellular solution containing 120 mM NaCl, 26 mM NaHCO 3 , 2.5 mM KCl, 1 mM NaH 2 PO 4 , 20 mM Glucose, 2.5 mM CaCl 2 , 1.3 mM MgSO 4 and adjusted to pH 7.4 and infused with 95% O 2 /5% CO 2 . Data was acquired using Igor Pro Software (Wavemetrics, Lake Oswego, OR). The stability of series and input resistances were confirmed throughout the experiment. Signals were filtered at 2 KHz and digitized at 5 KHz. To analyze action potential generation, cells were held in current clamp at 275 mV and a 10-500 ms current of 100-200pA was injected to depolarize cells to threshold.

Differentiation into Functional Neurons
Dental pulp cells were cultured and reprogrammed into iPSCs, then induced to differentiate into neurons, as described in the methods section. In previous experiments using neurons derived from F-iPSCs, the differentiation protocol resulted in the production of a heterogeneous mix of glutamatergic and GABAergic neurons that express forebrain, midbrain and hindbrain transcription factors (TFs) [1][2][3][4][5][6][7][8][9]. A small fraction of cells (,1%) express the dopamine marker TH (unpublished observations). TIPS4 produced a similar mix of neurons ( Figure 1A, 1B). As the cells matured, robust staining for the pre and post synaptic glutamatergic markers synaptophsin and PSD95 was seen ( Figure 1C). In addition, after exposure to a depolarizing current, a train of action potentials could be detected; sodium channel activation is responsible for the response, since it's blocked by tetrodotoxin (TTX) ( Figure 1D). These findings confirm that functional neurons can be developed from T-iPSCs.

RNA-Seq
Whole genome transcriptome analysis (RNA-Seq) was carried out on TIPS4 and TIPS4-C5 iPSCs, NPCs and 14 day neurons following neuronal differentiation. The expression levels for 1368 genes were differentially expressed during the transition from iPSCs to NPCs (627 higher in NPCs; 741 lower in NPCs, FDR ,0.05). A comparison between iPSCs and neurons showed 2543 differentially expressed genes (1286 increased in neurons; 1257 lower in neurons, FDR ,0.05) (entire list of differentially expressed genes during differentiation at all 3 transition points [iPSCs to NPCs; iPSCs to neurons; NPCs to neurons can be found in Table S2, Table S3, Table -S4). Among the genes that showed substantial decreases in expression in NPCs and neurons were POU5F1 (OCT4), LIN28A, TDGF1, which are expressed at high levels in embryonic and pluripotent stem cells, and are inactivated during differentiation (Table S2 and Table S3) [14,24].
Among the top coding genes that increased in expression during differentiation into NPCs and neurons were a number of members of the HOX gene family (discussed below), and the neuronal TFencoding genes POU3F3, MYT1L, and DLX1 ( Table 2, Table 3). MYT1l, along with POU3F2, and ASCL1, can reprogram fibroblasts directly to neurons [31]. POU3F2 and ASCL1 also increased in neurons, but just failed to meet the FDR criterion for genome-wide significance.
Several lncRNAs were also among the top genes that increased in expression. However, expression levels were relatively low (FPKM ,5), with the exception of LINC00473, a long intergenic non-coding RNA that overlaps with RP11-252P19.3 (see below).
Finally, there were 289 genes that showed a significant difference in expression during the transition from NPCs to neurons (168 increased; 121 decreased) (Table S4). Interestingly, several lncRNAs showed a significant increase in expression in neurons compared with NPCs, suggesting that they could be involved in neuronal maturation, rather than the initial differentiation of iPSCs into NPCs per se. These included RP11-466P24.7, RP11-64K12.10 and RP11-252P19. RP11-466P24.7 maps to the 39-UTR of an isoform of SV2C (synaptic vesicle glycoprotein 2C); RP11-64K12.10 maps near DISP2, which is involved in hedgehog signaling, and RP11-252P19.3 is embedded within SDIM1, which is down regulated in Alzheimer's brains and may affect NPC cell death [32,33].
Genes that increased or decreased in expression during differentiation were subjected to pathway analysis. As seen in  Table 4 and Table 5, enrichment for genes involved in SZ was found among those that increased in neurons, while genes involved in cell cycle regulation were decreased. A complete list of genes enriched for these and other functions can be found in Table S5 and Table S6. The findings show that iPSCs derived from dental pulp can be used for disease modeling neuropsychiatric disorders.

Comparison of T-iPSC and F-iPSC Neurons
RNA-seq profiles were obtained for neurons derived from the two T-iPSC lines and two F-iPSC lines. A total of 329 genes were differentially expressed at a nominally significant p-value (p,0.05), of which 63 remained significant after correcting for multiple testing (FDR ,0.05; 54 expressed at a lower level in the T-iPSCs; 9 expressed at higher levels, Table S7). qPCR was used to validate the RNA-seq findings for 8 differentially expressed genes in one tooth vs fibroblast set of neurons (ASCL1, EMX1, EMX2, FOXG1, LHX2, OTX2, TBR1, MYT1L, and FOXP2). The fold change differences were consistent with the RNA-seq findings (Table S8).
No significant differences were detected for any neurotransmitter receptor or transporter gene, with the exceptions of a significant decrease in the T-iPSC neurons in the level of GRID2 mRNA (glutamate receptor delta-2; a member of the ionotropic glutamate receptor family, and lower levels of SLC6A1 (vesicular GABA transporter), SLC5A7 (choline transporter) and SLC6A5 (glycine transporter) ( Table S7).
The most differentially expressed genes were SLITRK2 and SCUBE2. SLITRK2 codes for an integral membrane protein that shares homology with neurotrophin receptors; it has been implicated in a small subset of patients with BD and ASD [34,35]. SCUBE2 is expressed in the hindbrain and forms a complex with Sonic hedgehog and its receptor PTC1 to activate SHH-signaling [36,37].
The most striking difference between the tooth and skin-derived neurons is the substantially lower expression of a number of HOX genes (Table S7; Figure 2). As seen in the figure, which shows the relative expression of all HOX genes with FPKM .1, there were 14 that showed significantly lower levels in the tooth-derived neurons (FDR,0.05; depicted by an asterisk). Seven other HOX genes (HOXB-AS4, HOXA-AS3, HOXA5, HOXC8, HOXA9, HOXA7, and HOXA4) showed decreases that were nominally significant (p,0.05; FDR.0.05). HOX gene expression is involved in brain patterning and is regulated by retinoic acid (RA), among other signaling pathways. Interestingly, several others genes regulated by RA are expressed at significantly lower levels in the tooth samples, including RBP4 and RORB, and the homeobox genes PHOX2A, IRX1, and IRX2 (Table S7).
To examine differentially expressed genes more systematically, we subjected the data to Ingenuity Pathway Analysis (IPA). Interestingly, neurological disease/schizophrenia was the top category for genes that were expressed at higher levels in the neurons derived from teeth (nominally significant p-value of ,0.05) ( Table 6).

Discussion
Disease modeling using iPSCs must be carried out using readily accessible sources of somatic cells from patients for reprogramming, such as skin fibroblasts, hair keratinocytes, CD34+ leukocytes, epithelial cells found in urine, and dental pulp [11,[38][39][40][41]. Of these, a skin biopsy is the most invasive, which would make it the least suitable for children with developmental disorders. Obtaining hair follicles and blood for iPSC reprograming are rather non-invasive, certainly, but are not totally free of causing some degree of distress in autistic and developmentally disabled children. Deciduous teeth that are naturally shed during childhood, on the other hand, provide a source of cells for reprogramming that is not invasive, posing no additional stress to the child. Dental pulp, however, would appear to be the least convenient for researchers, since it relies on waiting for deciduous teeth to be shed. Yet, considering the time required to generate and characterize iPSC lines and the fact that children lose 20 deciduous teeth between the ages of ,5-12, with minimal planning, collecting a library of dental pulp cells for iPSC reprogramming should not be a limiting factor. From a biological perspective, dental pulp could prove to be a better source of iPSCs for disease modeling neuropsychiatric disorders because of its developmental origins. During early vertebrate development, embryonic ectoderm differentiates into neural and neural plate borders, and epidermal regions [42,43]. Dental pulp contains ectomesenchyme, which is derived from ectoderm, specifically neural crest cells, while fibroblasts are derived from ectoderm that is programmed to become epidermis [42,43]. Considering the fact that gene expression could be affected by the retention of some epigenetic marks following reprogramming that are dependent on the somatic cell of origin, an assessment of gene expression profiles using neurons derived from different reprogrammed cells could show differences that might be relevant to in vitro disease modeling, a question we have addressed in this paper. Expression profiling showed that neurons derived from T-iPSCs and F-iPSCs differed for some key genes, notably multiple members of the HOX gene families. HOX gene expression is involved in anterior/posterior patterning and the development of hindbrain structures. The homeobox genes IRX1 and IRX2, which are also involved in brain patterning, were expressed at significantly lower levels in the neurons from T-iPSCs as well. Since both the HOX and IRX gene families are induced by RA, lower levels of expression could reflect a lower sensitivity to the RA present in the medium used during the development of NPCs.
While lower levels of expression for genes involved in hindbrain development were seen in the neurons derived from teeth, several TFs involved in the forebrain development were significantly increased as well; most notably FOXP2, OTX1, and LHX2. FOXP2 codes for a TF involved in the development of communication and language neural networks that has been implicated in ASD [44,45]. Considering the fact that SZ and ASD are associated with cognitive abnormalities, the decrease in expression of genes involved in hindbrain development and an increase in expression of some key forebrain TFs suggests that neurons derived from T-iPSCs may have some advantages over those derived from fibroblasts in conditions like SZ and ASD that are associated with cognitive and language impairment. In fact, pathway analysis of differentially expressed genes showing enrichment for genes involved in SZ supports this notion. On the other hand, there may be a disadvantage for disorders affecting hindbrain structures. Whether the differences in gene expression persist using other neuronal differentiation methods remains to be seen. Figure S1 A. Immunocytochemistry for pluripotency markers (Tra-1-60, Tra-1-80, SSEA3, SSEA4) and DAPI nuclear stain (blue) for clone TIPS4. B. Expression of germ layer markers; AFP (endoderm), desmin (mesoderm) and b-IIItubulin (ectoderm). (TIF) Figure S2 Immunocytochemistry for pluripotency markers clone TIPS4-C5. A. Immunocytochemistry for pluripotency markers (Tra-1-60, Tra-1-80, SSEA3, SSEA4) and DAPI nuclear stain (blue) for clone TIPS4-C5. (TIF) Table S1 RNA-seq statistics. RNA-seq statistics for all T-IPSC and F-IPSC samples, and correlation coefficients for the T-IPSC samples (iPSCs, NPCs and neurons), and F-iPSCs. (XLSX) Table S2 Differentially expressed genes during transition from iPSCs to NPCs. Genes that significantly changed in expression during transition from iPSCs to NPCs for TIPS4 and TIPS4-C5 in descending order of significance. (XLSX)