FoxP1 orchestration of ASD-relevant signaling pathways in the striatum.

Mutations in the transcription factor Forkhead box p1 (FOXP1) are causative for neurodevelopmental disorders such as autism. However, the function of FOXP1 within the brain remains largely uncharacterized. Here, we identify the gene expression program regulated by FoxP1 in both human neural cells and patient-relevant heterozygous Foxp1 mouse brains. We demonstrate a role for FoxP1 in the transcriptional regulation of autism-related pathways as well as genes involved in neuronal activity. We show that Foxp1 regulates the excitability of striatal medium spiny neurons and that reduction of Foxp1 correlates with defects in ultrasonic vocalizations. Finally, we demonstrate that FoxP1 has an evolutionarily conserved role in regulating pathways involved in striatal neuron identity through gene expression studies in human neural progenitors with altered FOXP1 levels. These data support an integral role for FoxP1 in regulating signaling pathways vulnerable in autism and the specific regulation of striatal pathways important for vocal communication.

Autism spectrum disorder (ASD) denotes a group of heterogeneous neurodevelopmental conditions that are all characterized by diminished sociability, impaired communication, restricted interests, and stereotypic behaviors. While there is a strong genetic component to ASD, this is divided among several hundred genes, each with only a small contribution to the prevalence of the disorder (Geschwind and State 2015). Furthermore, many autism-risk genes are thought to exert their effects during early brain development (State and Sestan 2012; Xu et al. 2014; Parikshak et al. 2015). Transcription factors play a key role in orchestrating the spatial and temporal gene expression patterns important for this process. Therefore, the identification of gene networks regulated by transcription factors implicated in both ASD and brain development should provide insight into the complex developmental brain mechanisms at risk in autism.

The transcription factors Forkhead box P1 (FOXP1) and FOXP2 have been implicated in neurodevelopmental disorders such as ASD and developmental verbal dyspraxia (DVD), respectively (Bacon and Rappold 2012). Foxp1 is a member of the Fox family of transcription factors, for which there is a designated protein nomenclature (uppercase for primates, lowercase for rodents, and mixed case for other species) (Kaestner et al. 2000). Foxp1 is highly enriched within the developing and mature neocortex, hippocampus, and striatum (Ferland et al. 2003; Teramitsu et al. 2004). Numerous studies have identified heterozygous deletions, point mutations, and duplications of FOXP1 as being causal for ASD (Bacon and Rappold 2012). In particular, recent large-scale exome sequencing efforts have identified FOXP1 as a gene with recurrent de novo mutations associated with ASD (Iossifov et al. 2014). Therefore, understanding how FOXP1 functions within the brain should allow for key insights into the molecular pathways at risk in ASD. Several reports have begun to elucidate a role for Foxp1 in the brain (Rousso et al. 2012; Tang et al. 2012), and recent work has shown that mice with brain-specific loss of Foxp1 have altered hippocampal electrophysiology, striatal morphology, and social behaviors (Bacon et al. 2015). However, the region-specific transcriptional profile of Foxp1 in the mouse brain, how well this profile is conserved in human-relevant Foxp1 haploinsufficient models, and the behavioral consequences of disrupting these regional gene networks remain largely unknown.

FOXP2 is a paralog of FOXP1, and mutations in the FOXP2 gene lead to a number of brain and cognitive deficits, including DVD (Fisher and Scharff 2009; Bacon and Rappold 2012). In addition to being able to heterodimerize with Foxp2, Foxp1 expression overlaps with Foxp2 expression in the GABAergic medium spiny neurons (MSNs) of the striatum, a brain region critically involved in human language, vocal imitation in zebra finches, and rodent ultrasonic vocalizations (USVs) (Ferland et al. 2003; Li et al. 2004; Teramitsu et al. 2004; Fisher and Scharff 2009). Additionally, Foxp2 mutant mice demonstrate disruptions in mouse USVs as well as alterations in the electrophysiological and projection properties of MSNs (Shu et al. 2005; Enard et al. 2009; Vernes et al. 2011; French et al. 2012). Given the role for both Foxp1 and Foxp2 in the striatum, we hypothesized that Foxp1 regulates regional gene expression patterns in the brain and that normal levels of Foxp1 are crucial for mouse vocalization behavior. To test this hypothesis, we took advantage of a heterozygous (Foxp1+/−) mouse model and a human neural progenitor (hNP) cellular model with altered expression of FOXP1. Using high-throughput sequencing technologies, we used these two systems to identify a conserved role for FoxP1 in regulating autism-risk genes. We showed that Foxp1 differentially regulates the excitability of dopamine receptor 1-positive (D1+) versus D2+ MSNs. We also demonstrated reduced USVs in Foxp1+/− mice, similar to that seen in Foxp2+/− mice (Shu et al. 2005). This similarity in behavioral phenotype is reflected at the genomic level, as Foxp1-regulated genes in the striatum overlap with genes regulated by Foxp2 in the striatum. Finally, we found that FoxP1 regulates conserved pathways involved in striatal identity in both humans and mice. Taken together, these results suggest that FoxP1 plays a critical role in regulating striatal function and vocal communication, which, when disrupted, contributes to phenotypes characteristic of ASD.

Results

Foxp1 gene regulation within distinct brain regions

In order to assess the ASD-relevant role of Foxp1 within the brain, we took advantage of a Foxp1 animal model. As Foxp1 knockout mice are embryonic-lethal at embryonic day 14.5 (E14.5) due to a developmental heart defect (Wang et al. 2004) and as most patients with FOXP1 mutations are haploinsufficient, we carried out analyses on Foxp1 heterozygous (Foxp1+/−) mice (Hu et al. 2006). We tested the specificity of an antibody recognizing FoxP1 using hNPs with forced FOXP1 expression as well as whole brains from E13.5 Foxp1 knockout embryos. We identified expression of two Foxp1 isoforms (A and D), previously shown to be expressed in mouse brains (Wang et al. 2003), both of which were absent in brain tissue from knockout embryos (Supplemental Fig. 1A). Three brain regions relevant to ASD with substantial levels of Foxp1 expression are the striatum, hippocampus, and neocortex (Ferland et al. 2003; Maloney et al. 2013). We quantitatively determined an ∼50% reduction in total Foxp1 protein levels (isoforms A and D) in the Foxp1+/− hippocampus or striatum compared with control littermates (Fig. 1A,B). Interestingly, neocortical expression of Foxp1 (either total protein or isoform-specific expression) was not reduced to 50% in Foxp1+/− mice, suggesting a homeostatic up-regulation of Foxp1 in the neocortex of these animals (Fig. 1A,B; Supplemental Fig. 1B,C).

Figure 1.

Regulation of ASD genes by Foxp1 in the mouse brain. (A) Representative immunoblot displaying reduced Foxp1 protein levels in the hippocampus (HIP) and striatum (STR), but not the neocortex (CTX), of Foxp1+/− mice. Gapdh was used as a loading control. (B) Quantification of Foxp1 expression in adult Foxp1+/− mouse brains. Data are represented as means ± SEM. n= 4 mice per genotype for each region. (*) P = 0.033 (hippocampus); (*) P = 0.0163 (striatum), Student's t-test, compared with wild-type levels normalized to Gapdh. (C) Venn diagram showing overlaps between the differentially expressed genes (DEGs) in the mouse and ASD gene lists (144 genes between the hippocampus and striatum [P = 1.21 × 10−26], 116 genes between the hippocampus and ASD [P = 3.74 × 10−9], and 43 genes between the striatum and ASD [P = 0.002], hypergeometric test [P-values were adjusted using Benjamini-Hochberg FDR procedure]). (D) Confirmation of salient ASD-related gene targets in independent striatal samples from Foxp1+/− mice using quantitative RT–PCR (qRT–PCR). Data are represented as means ± SEM. n = 4 mice per genotype. With the exception of Dner, all qRT–PCR values displayed are significant at P < 0.05 (Student's t-test, compared with wild-type levels normalized to actin). (E) Visualization of a striatal-specific submodule (MsM18) that contains Dpp10 (dipeptidyl peptidase) as a major hub gene. (F) qRT–PCR confirmation of Dpp10 and Kcnd2 activation in Foxp1+/− mouse striatal samples. Data are represented as means ± SEM. n = 4 mice per genotype. All qRT–PCR values displayed are significant at P < 0.05 (Student's t-test, compared with wild-type levels normalized to actin).

We ascertained potential transcriptional targets of Foxp1 in vivo using RNA sequencing (RNA-seq) in the hippocampus or striatum of Foxp1−/+ mice and control littermates. To identify differentially expressed genes (DEGs), we filtered using a false discovery rate (FDR) of <0.05 and an absolute log fold change of ≥0.3 (Supplemental Table 1). As a control, RNA-seq was also conducted in the neocortex, and, not unexpectedly, we did not observe significant changes in gene expression in this brain region (data not shown).

Foxp1 regulation of ASD-associated pathways in the striatum and hippocampus

To characterize the identified Foxp1 targets with respect to ASD etiology, we compared the list of Foxp1 DEGs with the current list of annotated ASD genes in the SFARI database (667 genes) (http://www.sfarigene.org). We found that in vivo Foxp1-regulated genes significantly overlap with ASD genes in both the hippocampus and striatum (Fig. 1C). The SFARI database stratifies genes based on the strength of their association with ASD, and when we removed the genes in categories #5 and #6 (hypothesized support and not supported, respectively) and repeated our analyses, we obtained a similar result (17 genes [P = 0.057] for striatum and 39 genes [P = 0.0001] for hippocampus, hypergeometric tests) (data not shown). Using quantitative RT–PCR (qRT–PCR), we confirmed 11 of 12 selected targets from the overlap between the Foxp1+/− striatal data set and the ASD genes in independent samples (Fig. 1D).

Mouse Foxp1 targets were further prioritized with respect to neurodevelopmental human diseases using weighted gene coexpression network analysis (WGCNA), which allows for the discovery of networks (or modules) of genes with high levels of coexpression (Supplemental Fig. 2A; Zhang and Horvath 2005; Oldham et al. 2008). The top hub gene (or gene with the highest number of connections) in the striatal-associated MsM18 module is Dpp10 (dipeptidyl peptidase) (Fig. 1E). DPP10 is an ASD gene that encodes for a protein that regulates surface expression and properties of the potassium channel Kv4.2 (Marshall et al. 2008; Foeger et al. 2012). Of note, the gene encoding Kv4.2, KCND2, has also been implicated in ASD (Lee et al. 2014) and is highlighted within the MsM19 module (Supplemental Fig. 2C). We also observed and confirmed that Dpp10 is increased and that Kcnd2 is decreased in the striatum of Foxp1+/− mice (Fig. 1F).

Alteration of Kv4.2 function has been previously observed in a mouse model of Fragile X syndrome (FXS) (Gross et al. 2011), and Fragile X mental retardation protein (FMRP)-regulated genes have previously been shown to have significant genomic interactions with ASD-relevant pathways in human brain development (Parikshak et al. 2013). We therefore compared the mouse WGCNA modules with previously identified FMRP targets (Darnell et al. 2011) and found modules containing FMRP targets (MsM1, MsM6, MsM12, MsM14, and MsM23) (Supplemental Fig. 2A). While certainly interesting with regard to potential converging pathways, such enrichments need to be interpreted cautiously, as recent work has uncovered that FMRP targets tend to be highly expressed long genes in the brain (Ouwenga and Dougherty 2015). Of particular interest, MsM14 correlates with genotype and contains a number of FMRP target genes, including the gene encoding FMRP (Fmr1) (Supplemental Fig. 2B), highlighting a potential direct role for coordination of disease-relevant genes in the striatum by Foxp1 and FMRP.

Foxp1 regulates shared targets with Foxp2 in the striatum

As previous work has implicated a role for the related transcription factor FoxP2 in striatal function, including altered MSN electrophysiology and morphology (Enard et al. 2009), and as the striatum is also one of the few brain regions where Foxp1 and Foxp2 have overlapping expression (Ferland et al. 2003), we compared the list of Foxp1 target genes in the striatum with published striatal Foxp2 targets in Foxp2+/− mice (Enard et al. 2009). We identified a significant overlap between Foxp1-regulated genes and previously published Foxp2 targets that are changing in the same direction across data sets with reduction of the respective transcription factors, indicating possible coregulation of these targets (Fig. 2A). This overlap represents 12% of the total Foxp1 target genes identified in the striatum. Using independent samples, we confirmed six of these genes changing with Foxp1 expression in the striatum via qRT–PCR (Fig. 2B). Within the in vivo WGCNA analysis, both Foxp1 and Foxp2 are coexpressed within the MsM3 module, which is enriched for striatal DEGs (Fig. 2C). Interestingly, within the MsM3 module, the gene encoding the dopamine receptor Drd1a is coexpressed with both Foxp2 and Foxp1 (Fig. 2C).

Figure 2.

Foxp1 and Foxp2 regulate overlapping targets within the striatum. (A) Significant overlap of DEGs in the striatum of Foxp1+/− and Foxp2+/− mice (67 genes between the Foxp1+/− and the Foxp2+/− striatal data sets [P = 2.82 × 10−5], hypergeometric test). (B) qRT–PCR confirmation of a subset of these genes in independent Foxp1+/− striatal samples. Data are represented as means ± SEM. n = 3 mice per genotype. All qRT–PCR values displayed are significant at P < 0.05 (Student's t-test, compared with wild-type levels normalized to actin). (C) Visualization of the regionally specific striatal module MsM3 showing coexpression of both Foxp1 and Foxp2. Foxp1 and Foxp2 connections are highlighted in magenta. Genes in bold typeface indicate striatal DEGs, and boxed genes indicate Foxp1 and Foxp2 DEGs that overlap. (D) RNA-seq data from Foxp1+/− mice and microarray data from Foxp2+/− mice were overlapped with the most recently published list of known enriched transcripts within D1+ or D2+ MSNs (Maze et al. 2014). Genes from both Foxp1+/− and Foxp2+/− mice significantly overlapped with D1+ MSN-enriched genes (36 genes [P = 1.12 × 10−5] and 61 genes [P = 1.99 × 10−12], respectively, hypergeometric test). P-values for each overlap are shown within bar graphs.

MSNs of the striatum are categorized as either D1+ (expressing the Drd1a receptor) or D2+ (expressing the Drd2 receptor) projection neurons, and these two subpopulations of neurons are associated with opposing functions in the coordination of motor activity (Gerfen and Surmeier 2011). To investigate whether disrupted Foxp1 signaling in the striatum would be expected to produce differential gene expression changes in D1+ versus D2+ MSNs, we overlapped our RNA-seq data set with published gene lists obtained from translating ribosome affinity purification of D1+ and D2+ MSNs (Maze et al. 2014). We found a significant enrichment of both Foxp1 and Foxp2 target genes within D1+ MSNs specifically (Fig. 2D). Although we found that the number of Foxp1 target genes is roughly equally distributed between genes enriched in both D1+ and D2+ MSNs (Fig. 2D), the number of Foxp2 target genes enriched in D1+ MSNs is almost twice the number of Foxp2 target genes enriched in D2+ MSNs (Heiman et al. 2008; Vernes et al. 2011; Maze et al. 2014). These results are in line with published data showing that Foxp2 is more enriched in D1+ MSNs. Moreover, the overlapping targets of Foxp1 and Foxp2 going in the same direction are expressed only in D1+ MSNs, supporting coordinated regulation in these specific neurons. Interestingly, Foxp1-specific target genes that are enriched in D2+ MSNs include several genes involved in cation transport (e.g., Atp1b1, Kcnk2, Htr7, Kcnip2, and Hrh3) (Supplemental Table 2). Together, these data support a role for Foxp1 and Foxp2 providing coordinated regulation of striatal signaling pathways and that this regulation may be differential between D1+ and D2+ MSNs in Foxp1+/− mice.

Reduction of Foxp1 leads to differential changes in the excitability of striatal MSNs

Together with the coregulation of genes by Foxp1 and Foxp2 (Fig. 2), the gene expression data indicated a role for FoxP1 in regulating genes coding for proteins involved in both ion channel and neuronal activity, in particular within D2+ MSNs (Supplemental Tables 1, 2). We therefore investigated the effect of reduced Foxp1 expression on neuronal activity within either D1+ or D2+ MSNs. At postnatal day 18 (P18), acute striatal slices were made from progeny of Foxp1+/− mice crossed with either Drd1a-tdTomato+/− or Drd2-GFP+/− reporter mice (Gong et al. 2003; Ade et al. 2011) and whole-cell recordings of MSNs were carried out. D2+ (GFP+) MSNs from Drd2-GFP+/− mice (Fig. 3A) exhibited significantly increased excitability, as indicated by the higher number of action potentials evoked for a given current step (Fig. 3B,C), an increase in input resistance (Fig. 3D), and a decrease in current threshold (Fig. 3E). We also observed no differences in resting potential (Fig. 3F), the action potential width (Fig. 3G), or the frequency of the spontaneous excitatory postsynaptic events (sEPSCs) of D2+ MSNs (Fig. 3H). Interestingly, we did observe a significant decrease in the amplitude of sEPSCs of these MSNs (Fig. 3I). Together, these data demonstrate that reduction of Foxp1 leads to increased excitability of D2+ MSNs in response to reduced Foxp1 expression.

Figure 3.

D2+ MSNs of Foxp1+/− mice have increased excitability. (A) Example image of a recorded GFP+ (D2+) neuron. (B) Example recordings depicting spiking in response to a 125-pA current step in control and Foxp1+/− MSNs. (C) Firing rate versus input curves is significantly increased in Foxp1+/− MSNs. Data are represented as means ± SEM. n = 18 wild-type cells and 29 Foxp1+/− cells. (*) P = 0.040, two-way ANOVA with repeated measures for current step, compared between genotypes. (D) Input resistance is significantly increased Foxp1+/− MSNs. Data are represented as means ± SEM. n = 19 wild-type cells and 30 Foxp1+/− cells. (***) P = 0.0004, Student's t-test, compared between genotypes. (E) The minimum, threshold current required for evoking an action potential is significantly decreased in Foxp1+/− MSNs. Data are represented as means ± SEM. n = 19 wild-type cells and 30 Foxp1+/− cells. (*) P = 0.049, Student's t-test, compared between genotypes. (F) Resting potential is not significantly changed in Foxp1+/− MSNs. Data are represented as means ± SEM. n = 19 wild-type cells and 30 Foxp1+/− cells. P = 0.53, Student's t-test, compared between genotypes. (G) Action potential width is not significantly altered in Foxp1+/− MSNs. Data are represented as means ± SEM. n = 19 wild-type cells and 30 Foxp1+/− cells. P = 0.57, Student's t-test, compared between genotypes. (H) Spontaneous EPSC frequency is not significantly changed in Foxp1+/− MSNs. Data are represented as means ± SEM. n = 17 wild-type cells and 25 Foxp1+/− cells. P = 0.091, Student's t-test, compared between genotypes. (I) Spontaneous EPSC amplitude is significantly decreased in Foxp1+/− MSNs. Data are represented as means ± SEM. n = 17 wild-type cells and 25 Foxp1+/− cells. (**) P = 0.004, Student's t-test, compared between genotypes.

Given the opposing functions traditionally associated with D1+ and D2+ MSNs (Gerfen and Surmeier 2011) and the possibility for differential regulation of gene expression within D1+ and D2+ MSNs in the Foxp1+/− mouse striatum (Fig. 2D), we asked whether the increased excitability of D2+ neurons due to Foxp1 loss was generalizable to all MSNs. Again, at P18, we carried out whole-cell recordings on MSNs from acute striatal slices. Although trending toward a decrease in excitability, D1+ (tdTomato+) MSNs from Drd1a-tdTomato+/− mice (Supplemental Fig. 3A) exhibited no significant change in their excitability compared with controls (Supplemental Fig. 3B,C). We also found no significant increase in input resistance (Supplemental Fig. 3D) or current threshold (Supplemental Fig. 3E) and no significant difference in resting potential (Supplemental Fig. 3F) or action potential width with reduction of Foxp1 in these neurons (Supplemental Fig. 3G). Finally, we observed no changes in the frequency or amplitude (Supplemental Fig. 3H,I) of sEPSCs. These data indicate that haploinsufficiency of Foxp1 causes differential changes in the membrane excitability of D1+ and D2+ MSNs.

Foxp1 regulates mouse USVs

Huntington's and Parkinson's disease mouse models provide evidence for the involvement of MSNs in directing the production of USVs (Pietropaolo et al. 2011; Grant et al. 2014). Additionally, knockout of the Drd2 receptor reduces the number of USVs produced by mouse pups (Curry et al. 2013). Because we uncovered a significant overlap between Foxp1+/− and Foxp2+/− striatal target genes as well as altered MSN excitability as a response to loss of Foxp1, we hypothesized that reduction of Foxp1 would lead to an altered USV phenotype similar to that seen in Foxp2 mutant mice (Shu et al. 2005). To test this hypothesis, we examined USVs in a maternal separation paradigm. Paralleling what has previously been seen in Foxp2+/− mice (Shu et al. 2005), we observed a significant decrease in both the number of times a Foxp1+/− mouse pup called (“bouts”) (Fig. 4A) and the total number of calls (Fig. 4B) compared with littermate controls at P4 and P7 (see the Materials and Methods; Supplemental Fig. 4A for analysis details). Additionally, as a trend, the call bouts and total number of calls produced by the Foxp1+/− mouse pups are reduced across all days (Fig. 4A,B). We also observed a significant decrease in the mean call frequency, as a trend across all days, in the Foxp1+/− mouse pups (Fig. 4C). Other parameters, such as average call duration and the fraction of calls with jumps, were not altered (Fig. 4D,E). Interestingly, we observed that the average slope of a call was significantly decreased in Foxp1+/− mice compared with controls (Fig. 4F). This result is the opposite of the increase in call slope exhibited by humanized Foxp2 mice (Enard et al. 2009).

Figure 4.

Foxp1 haploinsufficiency results in reduced mouse vocalizations. (A) Foxp1+/− mouse pups exhibit a significantly reduced number of vocalization bouts. Data are represented as means ± SEM. n = 57 wild-type pups and 34 Foxp1+/− pups. (*) P = 0.033 at P4; (***) P = 0.0003 at P7, two-way ANOVA with a Sidak multiple comparison test, compared between genotypes. (B) Foxp1+/− mouse pups exhibit fewer total numbers of USVs at P7. Data are represented as means ± SEM. n = 57 wild-type pups and 34 Foxp1+/− pups. (*) P = 0.038 at P4; (**) P = 0.006 at P7, two-way ANOVA with a Sidak multiple comparison test, compared between genotypes. (C) As a trend, Foxp1+/− mouse pups exhibit a significant reduction in their mean call frequency across all days. Data are represented as means ± SEM. n = 57 wild-type pups and 34 Foxp1+/− pups. Two-way ANOVA with a Sidak multiple comparison test, compared between genotypes. (D) Foxp1+/− mice show no differences in average call duration. Data are represented as means ± SEM. n = 57 wild-type pups and 34 Foxp1+/− pups. P = 0.99, two-way ANOVA with a Sidak multiple comparison test, compared between genotypes. (E) Foxp1+/− mice show no difference in the fraction of calls with frequency jumps. Data are represented as means ± SEM. n = 57 wild-type pups and 34 Foxp1+/− pups. P = 0.27, two-way ANOVA with a Sidak multiple comparison test, compared between genotypes. (F) Foxp1+/− mice display a significant reduction in the average slope of a call at P10. Data are represented as means ± SEM. n = 57 wild-type pups and 34 Foxp1+/− pups. (**) P = 0.001, two-way ANOVA with a Sidak multiple comparison test, compared between genotypes. The main effects for genotype and postnatal day and the interactions between these two variables are reported at the bottom of each panel.

Differences in weight gain have been proposed to explain some variation seen in the postnatal USVs of transgenic mouse models (Scattoni et al. 2009). However, there were no significant differences in the weight gain of Foxp1+/− mice compared with controls (Supplemental Fig. 4B). To assess whether the vocalization deficits observed in the Foxp1+/− mice are secondary to a generalized impairment in striatal-mediated behaviors, we assessed locomotion in the open field test as well as rotorod performance, forelimb and hindlimb grip strength, nest building, and grooming behaviors in these animals (Supplemental Fig. 5). We also performed postnatal righting reflexes as part of an abbreviated SHIRPA battery to evaluate overall neurological function in these mice (see the Materials and Methods; Supplemental Figs. 5A, 6). In summary, we found that Foxp1+/− mice display no differences in either the SHIRPA test, righting reflexes, nest building, rotorod performance, or grooming behaviors. Interestingly, Foxp1+/− mice do display hyperactivity in the open field test and decreased performance in the forelimb and hindlimb grip test. Together, these data suggest that wild-type levels of Foxp1 expression are important for normal mouse vocal behavior but are not required for most striatal-based behaviors.

FOXP1 gene regulation in human neural cells

In order to identify FoxP1 target genes that are most relevant to human brain development and ASD, we characterized the FOXP1 target genes in hNPs, which demonstrate a higher fidelity with in vivo brain transcriptomic data than either human embryonic stem cells or induced pluripotent stem cells (Konopka et al. 2012; Stein et al. 2014) and are genetically tractable using lentiviruses (Konopka et al. 2009). As undifferentiated hNPs do not express FOXP1 endogenously and given the current paucity of chromatin immunoprecipitation (ChIP)-grade antibodies against FoxP1, we ascertained direct FOXP1 targets by transducing hNPs with a lentivirus containing Flag-tagged FOXP1 or a GFP control virus (Fig. 5A). Forced expression of FOXP1 was limited to the nucleus (Fig. 5D).

Figure 5.

Gene regulation by FOXP1 in human neural cells. (A) Representative immunoblot depicting overexpression of FOXP-Flag signal in hNPs transduced with a FOXP1-Flag expression construct (hNPFOXP1) but not in hNPs with a GFP expression construct (hNPGFP). β-Tubulin was used as a loading control. (B) Representative immunoblot confirming expression of FOXP1-Flag in input samples and enrichment of FOXP1-Flag during the immunoprecipitation (IP) portion of ChIP from hNPFOXP1 lysates. (C,D) Representative images of hNPGFP and hNPFOXP1 demonstrate that FOXP1 expression (red) in hNPFOXP1 is restricted to the nucleus (DAPI, blue) and that FOXP1 is not expressed within neurites (Tuj1, green) and is absent in hNPGFP. (E) Significant overlap between gene targets from RNA-seq and ChIP-seq (ChIP followed by DNA sequencing) performed on hNPFOXP1 (92 genes between hNPFOXP1 RNA-seq and hNPFOXP1 ChIP-seq [P = 4.43 × 10−5, hypergeometric test]). (F) Significant overlap among RNA-seq DEGs, ASD genes, and FMRP targets (102 genes between hNPFOXP1 RNA-seq and ASD genes [P = 0.013], 122 genes between hNPFOXP1 RNA-seq and FMRP genes [P = 0.023], and 125 genes between ASD and FMRP genes [P = 1.34 × 10−35], hypergeometric test [P-values were adjusted using Benjamini-Hochberg FDR procedure]). (G) qRT–PCR confirmation of a subset of these overlapping genes in independent hNPFOXP1 samples. Data are represented as means ± SEM. n = 4 samples per treatment. All qRT–PCR values displayed are significant at P < 0.05 (Student's t-test, compared with hNPGFP levels normalized to actin). (H) DEGs from these overlaps are equally represented among repressed and activated genes. (I, left panel) Human genome browser view showing the ChIP-seq result of enrichment of FOXP1 binding compared with GFP control. (Right panel) ChIP-PCR confirmation of enriched binding of DPP10 by FOXP1 in hNPFOXP1 compared with hNPGFP using two separate primer pairs (DPP10 primers A and B) compared with control primers. Quantified data are represented as means ± SEM, four samples per treatment. All qRT–PCR values displayed are significant at P < 0.05 (Student's t-test, compared with hNPGFP levels normalized to actin). (J) DPP10 is repressed with FOXP1 overexpression in hNPFOXP1 samples. Quantified data are represented as means ± SEM, four samples per treatment. All qRT–PCR values displayed are significant at P < 0.05 (Student's t-test, compared with hNPGFP levels normalized to actin).

To identify genome-wide direct targets of FOXP1, we conducted both RNA-seq and ChIP followed by DNA sequencing (ChIP-seq) in hNPs overexpressing FOXP1. Using the Flag tag on FOXP1, we identified >600 genes enriched for FOXP1 binding (Fig. 5B,E; Supplemental Table 1). These directly bound targets are enriched for forkhead motifs (Supplemental Fig. 7; Stroud et al. 2006). Again, using RNA-seq, an FDR of <0.05, and an absolute log fold of ≥0.3, we uncovered >1500 DEGs within this cellular paradigm (Fig. 5E; Supplemental Table 1). These DEGs are significantly enriched for gene ontology (GO) categories such as axon guidance, neuronal development, and neuronal differentiation, and the overlap between both ChIP and RNA-seq data represents directly regulated FOXP1 targets in hNPs (Supplemental Table 3). RNA-seq and ChIP-seq genes significantly overlap (Fig. 5E); however, because this overlap is significant yet small, these results suggest that the majority of gene regulation by FOXP1 occurs through indirect effects on signaling cascades, as might be expected for a transcription factor.

FOXP1 regulates ASD-relevant genes in hNPs

To further characterize the identified hNP FOXP1 targets with respect to ASD etiology, we again compared the list of FOXP1 DEGs to the current list of annotated ASD genes in the SFARI database. We observed a significant overlap of FOXP1 targets and ASD genes (Fig. 5F). When we overlapped the list of hNP DEGs with the curated list of ASD genes (i.e., not including genes in categories #5 or #6), we also obtained a significant overlap (48 genes, P = 0.023, hypergeometric test) (data not shown). hNP DEGs that overlapped with the SFARI gene database were selected and confirmed within an independent hNP cell line using an independent measure of expression: qRT–PCR (Fig. 5G). Previous work suggested that the members of the Foxp subfamily of forkhead transcription factors are primarily transcriptional repressors (Wang et al. 2003). However, we showed that the related transcription factor FOXP2 is also able to activate transcription (Spiteri et al. 2007). In line with those data, we found an almost equal representation of activated and repressed FOXP1 targets that overlap with ChIP-seq and ASD lists (Fig. 5H). Additionally, we also confirmed that FOXP1 directly binds within the first intron of DPP10 and represses its expression in hNPs overexpressing FOXP1 (Fig. 5I,J). Together with the results from the Foxp1+/− mice (Fig. 1F), this indicates that Dpp10 is a conserved direct repressed target of FoxP1. Moreover, many genes overlapped with directional consistency between striatum and hNPs (12%) (Supplemental Table 1).

Using WGCNA again, we uncovered nine modules with first principal components correlating to FOXP1 expression (hNPM2, M3, M4, M6, M7, M13, M16, M20, and M21) (Supplemental Fig. 8). We then compared the hNPFOXP1 RNA-seq data to recently reported coexpression modules derived from in vivo developing human brains (Parikshak et al. 2013). We found a significant overlap of DEGs in the hNPFOXP1 data set with this report's M17 module (Supplemental Table 1). The M17 module is one of three modules previously identified to contain a significant overlap with known ASD genes. We also compared the hNPFOXP1 RNA-seq data with two other coexpression modules: asdM12 and asdM16 (derived from human brain tissue samples from ASD cases and controls), which were highly correlated with ASD disease status (Supplemental Table 1; Voineagu et al. 2011). We found that many genes within these two modules were also found within the modules correlating to FOXP1 expression. Interestingly, DPP10 is also present in asdM12, which further emphasizes its relevance to ASD etiology. Thus, the data from manipulation of FOXP1 expression in the in vitro system recapitulate identified genomic relationships from in vivo human brain data.

Conserved regulation of FoxP1 targets within the striatum

To further demonstrate the relevance of the Foxp1+/− mouse data with human biology and disease, we performed an analysis of module preservation (Langfelder et al. 2011) between either the Foxp1+/− mouse hippocampal or striatal WGCNA data and the hNP WGCNA data. This approach allows one to determine how conserved gene coexpression relationships are between the two species. Interestingly, we found that there was significantly greater preservation of modules between the Foxp1+/− mouse striatal modules and the hNP modules than between the Foxp1+/− mouse hippocampal modules and the hNP modules (Fig. 6A). To examine whether any of the preserved human coexpression modules contain specific transcription factor-binding motifs, we used the ChIP enrichment analysis (ChEA) database, which contains experimental ChIP and ENCODE data sets (Lachmann et al. 2010). We found enrichment of FOXP2 motifs as well as other autism-related transcription factors (Fig. 6B). Finally, we used a recently developed tool—cell-specific expression analysis (Xu et al. 2014)—to examine within which brain regions and cellular populations the conserved FoxP1 targets are enriched. We found that DEGs down-regulated with loss of Foxp1 and up-regulated with overexpression of FOXP1 are enriched for striatal genes (Fig. 6C). In contrast, genes up-regulated with loss of Foxp1 and down-regulated with overexpression of FOXP1 are enriched for neocortical genes (Fig. 6D). Together, these data suggest that FoxP1 regulates conserved pathways in both humans and mice that are important in preserving MSN identity.

Figure 6.

Coexpression network preservation between mouse and human data sets. (A) Module preservation analysis revealed that significantly more hNP modules are preserved in the striatum compared with the hippocampus. Zsummary scores >4 are well preserved, and those <2 are poorly preserved. (B) Genes in modules shared between humans and mice contain conserved binding sites for ASD-associated transcription factors, including FoxP2. (C) Genes down-regulated by loss of Foxp1 in mice and up-regulated by overexpression of FOXP1 in hNPs are enriched for striatal-associated genes. (D) Genes up-regulated by loss of Foxp1 in mice and down-regulated by overexpression of FOXP1 in hNPs are enriched for cortical genes. Briefly, hexagons are scaled to the stringency values of the specificity index thresholds (pSI), which ranks the region-specific enriched transcript gene lists from least specific to highly specific transcripts; i.e., outer hexagons represent larger, less specific lists (pSI of 0.05), while inner hexagons represent shorter, highly specific lists (pSI of 0.001). Bonferroni-Hochberg (BH)-corrected P-values are shown.

Discussion

Using unbiased genome-wide approaches in a patient-relevant Foxp1+/− mouse model and human neural cells, we uncovered a role for FoxP1 regulation of ASD-relevant genes. We observed that Foxp1 regulates gene expression in a region-specific manner within the brain, with the hippocampus and the striatum of Foxp1+/− mice containing DEGs enriched for distinct ontological categories. We also uncovered altered neuronal excitability in distinct populations of MSNs as well as gross alterations in the postnatal USVs of Foxp1+/− mice. Last, we provide evidence that FoxP1 regulates evolutionarily conserved neuronal pathways within the striatum, which are important for striatal identity.

The inclusion of FMRP target genes within FoxP1-correlated modules suggests overarching brain mechanisms at risk in ASD pathophysiology. FMRP is an RNA-binding protein that is expressed throughout the brain and is involved in dendritic morphology and plasticity through the translational regulation of numerous genes that function at the synapse (Darnell and Klann 2013). Deletions and mutations of the FMR1 gene can lead to FXS, which is characterized by autistic traits and intellectual disability (Hernandez et al. 2009). We found an enrichment of genes encoding ion channels altered in both the human and rodent FoxP1 models. For example, DPP10 is an ASD gene that is a conserved FoxP1 target (Figs. 1F, 5I,J). DPP10 functions to traffic surface expression of the KCND2 and KCND3 (or Kv4.2 and Kv4.3, respectively) potassium channels in neurons. We also uncovered activation of Kcnd2 by FoxP1 (Fig. 1F). Moreover, KCND2 is also an FMRP target (Kim et al. 2005); rare variants and genetic association of KCND2 have been reported in autism (Klassen et al. 2011; Lee et al. 2014), and impaired KCND2 function has been implicated in FXS (Gross et al. 2011). This convergence of Foxp1 downstream genes with FMRP-related genes suggests potential converging transcriptional and translational dysregulation in these disorders.

Relative to the entire brain, Foxp1 is among the top 100 enriched genes in the striatum (Heiman et al. 2008). This striatal enrichment of Foxp1 in the brain is greater than the comparative relative striatal expression of Foxp2. We showed significant overlaps between Foxp1 and Foxp2 gene targets in the striatum (in particular, D1+ MSN enriched genes) (Fig. 2A,D) and increased D2+ MSN excitability in Foxp1+/− mice (Fig. 3). Given that Foxp2 is preferentially expressed in D1+ MSNs (Vernes et al. 2011), we hypothesize that Foxp1 and Foxp2 may work in concert to differentially regulate neuronal excitability in these two populations of MSNs and therefore control striatal-based vocalizations. Therefore, the lack of alteration of D1+ MSN excitability in Foxp1+/− mice might be due to compensation by Foxp2, as supported by the overlapping target genes of Foxp1 and Foxp2 among D1+ MSN enriched genes. The identification of Foxp1-specific targets that are known to be involved in neuronal excitability within D2+ MSNs also supports the idea that differential gene regulation by Foxp1 in specific MSN subpopulations governs the observed neuronal and organismal phenotypes. For instance, Foxp1 may operate as a master regulator of genes important for overall neuronal function and activity in the striatum, with Foxp2 acting as a limiting factor for shared targets involved in vocalizations. This idea is bolstered by previous findings that Fmr1 and Foxp2 mutant mice exhibit increased striatal GABAergic transmission from and increased long-term depression in MSNs as well as decreased striatal volumes and deficits in postnatal USVs (Shu et al. 2005; Centonze et al. 2008; Enard et al. 2009; Roy et al. 2012; Ellegood et al. 2015). In addition, Foxp2 levels are unchanged in the striatum of Foxp1+/− mice (Supplemental Table 1) and do not appear to be significantly altered in either MSN population specifically (data not shown), further suggesting that alterations in Foxp2/Foxp1 stoichiometry in D1+ MSNs could be driving our findings. Finally, the significant gene coexpression module preservation between the mouse striatal and hNP gene expression data supports the relevance of these mouse data to a human disorder such as ASD. Given the evolutionary distance between these two species and the developmental differences between hNPs and the adult mouse striatum, it is remarkable that these correlations were found. Therefore, such a finding is evidence for the robustness and relevance of these gene coexpression networks with respect to FoxP1 expression and function.

These data suggest a role for Foxp1 in regulating ASD risk genes in a region-specific manner within the brain. In particular, we demonstrate that Foxp1 plays an important role in regulating genes involved in striatal development and function. In this study, we also provide the first evidence that Foxp1 specifically contributes to vocal communication. It will be important to determine how these changes occur throughout development in further experiments. Since the Foxp1+/− mice used in this study were whole-body knockouts and because Foxp1 has been shown to regulate the development of a host of organ systems (Wang et al. 2004; Hu et al. 2006; Shu et al. 2007; Dasen et al. 2008; Rousso et al. 2008), it cannot be entirely ruled out that the behavioral phenotypes displayed by these mice are secondary to the peripheral consequences of the knockout. Additionally, it should be noted that other brain regions besides the striatum, such as the neocortex, are known to both express Foxp1 and contribute to the production of USVs (Hisaoka et al. 2010; Sia et al. 2013). Moreover, while this study focused on a patient-relevant model of FOXP1 function (namely, haploinsufficiency), at least one study has demonstrated increased FOXP1 expression in lymphoblastoid cell lines from ASD patients (Chien et al. 2013). Therefore, the regional contribution and dosage relevance of FoxP1 to the behavioral manifestations presented in this study remain to be determined.

Materials and methods

Mice

Foxp1 heterozygous knockout (Foxp1+/−) mice were backcrossed with C57BL/6J mice for at least 10 generations to obtain congenic animals. Drd1a-tdTomato line 6 and Drd2-GFP reporter mice were generously provided by Dr. Craig Powell and maintained on a C57BL/6J background. Mice were kept in the barrier facilities of the University of Texas Southwestern Medical Center under a 12 h light–dark cycle and given ad libitum access to water and food. All studies with mice were approved by the University of Texas Southwestern Institutional Animal Care and Use Committee.

hNP cultures

hNP cultures were purchased from Lonza and maintained as previously described (Konopka et al. 2012). hNPs were transduced with lentiviruses containing pLUGIP-FOXP1-3XFlag or pLUGIP-GFP (control) and harvested 3 d after transduction for downstream applications, including immunoblotting, qRT–PCR, RNA-seq, and ChIP-seq.

RNA harvesting and real-time RT–PCR

RNA was purified from either hNPs or tissues dissected out from P47 male Foxp1+/− mice and littermate controls using an mRNeasy minikit (Qiagen) following the manufacturer's recommendations. qRT–PCR was performed as previously described (Spiteri et al. 2007). All primer sequences are available on request.

RNA-seq

mRNA was isolated from total RNA samples using polyA selection. Four independent samples from each brain region or cell type per genotype were included for a total of 24 mouse samples and eight human samples. Samples were randomized, and barcoded libraries were generated following the manufacturer's instructions (Illumina). RNA-seq was performed by the McDermott Sequencing Core at the University of Texas Southwestern Medical Center on an Illumina HiSeq 2000 sequencer (lllumina). Stranded, single-end 50-base-pair (bp) reads were generated for the hNP data, and stranded, paired-end 100-bp reads were generated for the mouse data.

RNA-seq data analysis

Reads were aligned to either hg19 or mm10 using TopHat (Trapnell et al. 2009) and Bowtie (Langmead et al. 2009). To obtain the gene counts, we used the HTSeq package (Anders et al. 2014), and the reads were normalized using the RPKM (reads per kilobase per million mapped reads) method (Mortazavi et al. 2008) implemented in the RSeQC package (Wang et al. 2012). For further analysis, we performed a sample-specific RPKM filtering considering genes with RPKM values of 0.5 in treatments or controls. EdgeR (Robinson et al. 2010) was used to detect the DEGs in each species. We applied a filter of FDR of <0.05 and absolute log fold change of >0.3 for both the human and mouse data sets. We then reconstructed the human and mouse coexpression networks using the R package WGCNA (Langfelder and Horvath 2008). Modules were characterized using the biweight midcorrelation followed by signed network topology for both human and mouse data. Modules containing ≥30 genes were included in our analyses. For module visualization, we used the publically available VisANT software (Hu et al. 2013). To determine the reliability of the WGCNA module characterization and the DEGs, we performed a permutation test randomizing 1000 times the expression data associated with each gene, calculated the DEGs, and then applied the same module characterization. None of the permuted data showed similar module detection or different expression profiles compared with the observed data. We then considered the detected modules, the detected DEGs, and the further gene overlaps significantly different from random expectation (permutation test, P = 0.001). To infer the significance of the potential overlaps, we adapted a hypergeometric test. The resultant P-values were adjusted using the Benjamini-Hochberg FDR method (Benjamini and Hochberg 1995).

Foxp2 microarray analysis

Data from project GSE13588 (Enard et al. 2009) were downloaded from Gene Expression Omnibus (GEO) (http://www.ncbi.nlm.nih.gov/geo). Only Foxp2 heterozygous and matching control samples were selected for further analysis. Microarrays were analyzed using the R programming language and Bioconductor packages. We determined gene expression levels (robust multichip average [RMA] values) and MAS5 detection P-values from the probes using the “affy” library (Gautier et al. 2004). We considered the probe sets detected in at least one sample for a P < 0.05. Differentially expressed probe sets were then determined adapting the f-test function implemented in the “multtest” library (Pollard et al. 2005). The resulting P-values were then adjusted with the Benjamini-Hochberg method. Probe sets were considered differentially expressed for an adjusted P < 0.05.

GO analysis

GO analysis was carried out using DAVID (http://david.abcc.ncifcrf.gov). A category containing at least three genes and a corrected P-value of <0.05 (Benjamini-Hochberg method) was considered significant.

Antibodies

The following antibodies were used for immunoblotting (IB), immunoprecipitation (IP), or immunocytochemistry (ICC): anti-β-tubulin (rabbit, 1:10,000; abcam, 6046 [IB]), anti-Flag (mouse, 1:10,000 [IB/ICC], 10 µg [IP]; Sigma, F1804), anti-Foxp1 (rabbit, 1:5000 [IB], 1:1000 [ICC]) (Spiteri et al. 2007), anti-Gapdh (mouse, 1:5000; Millipore [IB]), and anti-Tuj1 (mouse, 1:1000 [ICC]; Covance, MMS-435P).

ChIP-seq

Fifty-million hNPs were used per experimental condition. Cells were fixed in 1% methanol-free formaldehyde for 10 min at room temperature and then quenched with glycine (125 mM final). Cells were washed twice in 1× cold PBS, resuspended in 10 mL of lysis buffer (50 mM HEPES-KOH at pH 7.5, 140 mM NaCl, 1 mM EDTA, 10% glycerol, 0.5% IGEPAL-CA630, 0.25% TritonX-100, 10 µL/mL protease inhibitor [PI] cocktail [Sigma], 7 µL/mL PMSF), and incubated for 10 min on ice. Pelleted cell nuclei were then resuspended in 1 mL of nucleus lysis buffer (200 mM NaCl, 1 mM EDTA, 0.5 mM EGTA, 10 mM Tris-HCl at pH 8.0, 10 µL/mL PI, 7 µL/mL PMSF) and incubated for 10 min on ice. Samples were sonicated in 300 µL of shearing buffer (1 mM EDTA, 0.5 mM EGTA, 10 mM Tris-HCl at pH 8.0, 0.1% SDS, 10 µL/mL PI, 7 µL/mL PMSF) using a Bioruptor (Diagenode) at 3-min intervals for a total of 12 min. Ten percent of volume from each sample was collected for input controls. One-hundred micrograms of precleared sheared chromatin and 1 µg of msFlag antibody were incubated overnight at 4°C while rotating. Magnetic IgG Dynabeads (Invitrogen) were washed three times with 5 mg/mL BSA solution in PBS and then incubated with sheared chromatin/antibody solution for 2 h at 4°C. Magnets were applied to samples at 4°C, and beads were washed with 500 µL of each of the following solutions supplemented with PI and rotated for 5 min at 4°C followed by magnetic separation: (1) low-salt wash buffer (0.1% SDS, 1% TritonX-100, 2 mM EDTA, 20 mM Tris-HCl at pH 8.0, 150 mM NaCl), (2) high-salt wash buffer (0.1% SDS, 1% TritonX-100, 2 mM EDTA, 20 mM Tris-HCl at pH 8.0, 500 mM NaCl), (3) LiCl wash buffer (0.25 M LiCl, 1% IGEPAL-CA630, 1% deoxycholic acid, 1 mM EDTA, 10 mM Tris-HCl at pH 8.0), and (4) TE buffer. After washes, beads were resuspended in elution buffer (50 mM Tris-HCl, 10 mM EDTA, 1% SDS) and incubated for 15 min at 65°C with vortexing every 2 min. Beads were magnetically separated, supernatant was collected, and cross-linking of all samples and inputs was reversed overnight at 65°C. DNA was purified using Qiagen MinElute columns and quantified using a Qubit Fluorometer. Sequencing was performed by the University of Texas Southwestern Medical Center McDermott Sequencing Core.

ChIP-seq data analysis

Reads were mapped to the human genome (hg19) using TopHat (Trapnell et al. 2009) and Bowtie (Langmead et al. 2009). The aligned reads were subsequently downsampled according to the lowest number of reads detected, whereas the potential duplicated reads were removed using the Picard package (http://broadinstitute.github.io/picard). The uniquely mapped reads were then analyzed using MACS (Zhang et al. 2008) for the detection of potential peaks. PeakSplitter (Salmon-Divon et al. 2010) was used to subdivide the larger peaks into smaller, more precise peaks using a height filtering of 0.7. The FOXP1 peaks were further compared with the GFP peaks applying a tag density ratio (TDR). For further analysis, we considered FOXP1 peaks with a TDR >2.0. The uncovered peaks were then annotated using the annotatePeaks function implemented in the HOMER package (Heinz et al. 2010).

Immunoblotting

Cellular lysates were obtained using lysis buffer containing 0.5% Nonidet P-40,1 mM PMSF, 0.1 mM Na3VO4, 50 mM NaF, 1 uM DTT, 2 µg/mL pepstatin, and 1 µg/mL leupeptin. Tissue samples were lysed in buffer containing 1% Igepal, 1 mM PMSF, 0.1 mM Na3VO4, 2 µg/mL pepstatin, and 1 µg/mL leupeptin. Protein concentrations were determined using a Bradford assay (Bio-Rad). A total of 35–45 µg of each sample was run and processed following standard protocols for both HRP-conjugated and fluorescent secondary antibodies.

Immunocytochemistry

hNPs were grown on glass coverslips and fixed with 4% PFA in PBS for 15 min and then washed with TBS at room temperature. Cells were permeabilized with TBS-T (0.4% Triton-X) for 15 min at room temperature and then washed with TBS at room temperature. Cells were treated with a blocking solution made of 3% normal donkey serum in TBS-T (0.2% Triton-X) for 30 min at room temperature. Cells were then incubated with primary antibodies diluted in blocking solution overnight at 4°C. Afterward, cells were rinsed with TBS, treated with secondary antibodies diluted in blocking solution for 1 h, and then rinsed with TBS, all at room temperature. Slides were imaged using a Zeiss Observer.Z1 inverted microscope and ZEN 2011 software.

Electrophysiology methods

Electrophysiology

Acute brain slices were prepared from Foxp1+/− and Foxp1+/+ mice crossed with either Drd1a-tdTomato or Drd2-GFP reporter mice (P17–P20) with the following procedure. Mice were anesthetized with 125 mg/kg ketamine and 25 mg/kg xylazine, and the brains were removed. Thalamocortical slices (Agmon and Connors 1991) 300 µm thick were cut at ∼4°C in dissection buffer, placed in ACSF for 30 min at 35°C, and slowly cooled over the next 30 min to 21°C. Whole-cell recordings were performed in the dorsal striatum, and cells were targeted with IR-DIC optics in an Olympus FV300 confocal microscope. Recordings were performed at 21°C. Data were collected with a 10-kHz sampling rate and a 3-kHz Bessel filter. Striatal neurons were identified by GFP or tdTomato fluorescence using confocal microscopy.

Electrophysiology solutions

ACSF contained 126 mM NaCl, 3 mM KCl, 1.25 mM NaH2PO4, 2 mM MgSO4, 26 mM NaHCO3, 25 mM dextrose, and 2 mM CaCl2. All slices were prepared in the following dissection buffer: 75 mM sucrose, 87 mM NaCl, 3 mM KCl, 1.25 mM NaH2PO4, 7 mM MgSO4, 26 mM NaHCO3, 20 mM dextrose, 0.5 mM CaCl2, and 1 mM kynurenate. All solutions were pH 7.4. ACSF was saturated with 95% O2/5% CO2. Unless stated otherwise, the pipette solution consisted of 130 mM K-Gluconate, 6 mM KCl, 3 mM NaCl, 10 mM HEPES, 0.2 mM EGTA, 4 mM ATP-Mg, 0.3 mM GTP-Tris, 14 mM phosphocreatine-Tris, and 10 mM sucrose. This was adjusted to pH 7.25 and 290 mOsm. The junction potential was ∼10 mV and was not corrected.

Ultrasonic vocalization recordings

Acquisition and processing

USVs were recorded from pups isolated from their dams at P4, P7, and P10. Pups were placed into clean plastic containers inside soundproof styrofoam boxes and recorded for 3 min. Recordings were acquired using an UltraSoundGate condenser microphone (Avisoft Bioacoustics, CM16) positioned at a fixed height of 20 cm above the pups and were amplified and digitized (∼20 dB gain, sampled at 16 bits, 250 kHz) using UltraSoundGate 416H hardware and Avisoft RECORDER software (Avisoft Bioacoustics). Sound spectrograms were prepared in MATLAB (50% overlapping, 512-point Hamming windows), resulting in 1.024-msec temporal resolution and 488.3-Hz spectral resolution. Spectrograms were band-pass filtered to 20–120 kHz and filtered for white noise. Positions of ultrasonic calls were determined automatically using a previously published method (Holy and Guo 2005).

Spectral and temporal measurements

Vocalization behavior occured in spurts of activity (“bouts”) separated by longer pauses. To quantify bouts of vocalization, spectrograms were segmented using a pause length of ≥0.25 sec, which was chosen based on the empirical distribution of pause times between calls. All intercall pauses <0.25 sec represent constituents of the same bout of vocalization. The means of the dominant frequency (“mean frequency”) as well as the duration time of individual calls were averaged over all calls by animal. The presence of instantaneous pitch jumps in calls was determined by a previously published method (Holy and Guo 2005), and the fraction of all calls containing such jumps was determined for each animal. The trend slope (in hertz per millisecond) of calls lacking instantaneous pitch jumps was determined by linear regression, and slopes were averaged over all calls by animal.

Statistics

Differences between genotypes on all measured features of vocalization were assessed using two-way analysis of variance, testing for main effects of genotype, day, and interaction of genotype by day. Post-hoc multiple comparisons were assessed using Sidak's procedure. Features of vocalization were considered independently.

Postnatal righting reflexes

Righting reflexes were assessed in P4, P7, and P10 Foxp1+/− and littermate control pups. In brief, pups were placed in a supine position on a clean, unobstructed surface, and the time taken to right onto all fours was measured using a stopwatch. A pup failed the test if its time to right exceeded 1 min. In such cases, the time was scored as 60 sec. Each pup received one trial at each postnatal time point.

Open field test

The open field assay was performed on adult Foxp1+/− and littermate control mice by individually placing each animal in a 16-in × 16-in Plexiglass box and allowing them to explore the arena for 5 min. Videos of each mouse were obtained and scored for average velocity of movement and total distance moved using the EthoVision XT software package (Noldus).

Rotorod test

Adult mice were placed on a textured drum within individual lanes of a Series 8 IITC Life Science rotorod. The drum was programed to accelerate from 4 to 40 rpm within a maximum time frame of 300 sec. Each mouse was positioned forward on the drum, and sensors detected the latency to fall, maximum revolutions per minute at fall, and total distance travelled for each mouse. Sensors were manually activated whenever a mouse made a full rotation holding onto the drum. Mice were tested for three consecutive days with four trials per day, separated by 20-min intervals.

Grip strength test

Forelimb and hindlimb grip strength was measured on adult mice using Chatillon Force Measurement equipment. Forelimbs or hindlimbs of each mouse were placed on a mesh wire meter and pulled away from the wire at constant force. Five consecutive measurements were recorded for both hindlimbs and forelimbs and averaged for a final grip strength measurement for each mouse.

Nesting behavior

Mouse nesting behavior was analyzed using a previously described approach (Deacon 2006). Briefly, adult mice were singly housed overnight with 3 g of intact nestlet material in a clean cage. After 16–18 h, the amount of unused nestlet material was weighed, and the nests formed were assessed to generate a nest quality score of 1–5 for each mouse.

Grooming behavior

Grooming behavior was assessed in adult mice by individually placing each mouse in a clean cage without nesting material and allowing them to habituate for 10 min. Afterward, grooming behaviors were recorded using an HDR-CX535 Handycam video camera (Sony), and videos were then manually scored based on the number of grooming bouts and total time spent grooming for 10 min.

SHIRPA

A modified SHIRPA behavioral screen from Rogers et al. (1997) was performed on adult mice. First, mice were individually placed in a viewing jar for 5 min. During this time, mice were scored for (1) body position (inactive [0], active [1], or excessive activity [2]), (2) tremors (absent [0] or present [1]), (3) palpebral closure (open [0] or closed [1]), (4) coat appearance (tidy and well-groomed coat [0] or irregularities/piloerection [1]), (5) skin color (blanched [0], pink [1], or deep red [2]), (6) whiskers (absent [1] or present [0]), (7) lacrimation (absent [0] or present [1]), (8) defecation (absent [0] or present [1]), (9) gait (fluid with 3-mm pelvic elevation [0] or lack of fluidity [1]), (10) tail elevation (dragging [0], horizontal elevation [1], or elevated tail [2]), and (11) startle response (none [0], Preyer reflex [1], or reaction in addition to Preyer reflex [2]). Mice were then transferred to a clean cage, and the following behaviors were recorded in or above this arena: (12) touch escape (no response [0], response to touch [1], or flees prior to touch [2]), (13) trunk curl (absent [0] or present [1]), (14) limb grasping (absent [0] or present [1]), (15) pinna reflex (absent [0] or present [1]), (16) corneal reflex (absent [0] or present [1]), (17) contact righting reflex (absent [0] or present [1]), (18) evidence of biting (none [0] or biting in response to handling [1]), (19) vocalizations (nonvocal [0] or audible in response to handling [1]), (20) positional passivity (struggles when held by tail [0], when held by neck [1], or laid supine [2] or no struggle [3]). Both pinna and corneal reflexes were tested with a 0.15-mm-diameter nylon filament from Touch Test Sensory Evaluators (Semmes-Weinstein Monofilaments).

Other statistics

P-values were calculated with Student's t-test (two-tailed, type 2). F-values were calculated with two-way ANOVA followed by a Tukey post-hoc test for multiway comparison. Data were assumed to be normally distributed. P-values for overlaps were calculated with a hypergeometric test using a custom-made R script. We obtained an independent background for population size (for humans, human protein-coding genes [20,389 genes] and BrainSpan-expressed genes [15,585 genes] [Kang et al. 2011], and for mice, Allen brain-expressed genes [13,600 genes] [Lein et al. 2007]). We used the protein-coding genes for background in the hypergeometric test used in Figure 5E. We used the BrainSpan-expressed genes for background in the hypergeometric test used in Figure 5F and Supplemental Figure 8. We used the Allen brain-expressed genes for background in the hypergeometric test used for Figure 1C and 2, A and D, and Supplemental Figure 2A. P-values were adjusted for multiple comparisons using the Benjamini-Hochberg FDR procedure when required. A two-way permutation test of 1000 was adapted to validate the overlaps. First, we randomized the external gene sets (for example, ASD or FMRP) by randomly selecting the same number of genes from an independent brain-expressed gene list (for humans, BrainSpan-expressed gene list; for mice Allen-expressed gene list) and subsequently calculating the overlap P-values. The second approach randomized the internal gene sets (for example, STR_DEG or hNP_DEG) by randomly selecting the same number of genes from RNA-seq-expressed genes and subsequently calculating the overlap P-values. Moreover, we adapted a permutation test to evaluate the detected DEGs, randomizing 1000 times the RNA-seq data and recalculating the DEGs. Analysis for RNA-seq, ChIP-seq, and microarrays were performed using custom-made R scripts implementing functions and adapting statistical designs comprised in the libraries used.

Accession Numbers

The NCBI GEO accession number for the next-generation sequencing data reported in this study is GSE62718.