Highly conserved upstream sequences for transcription factor genes and implications for the regulatory network

Results

Alignment of the 8-kb Upstream Sequences of the Human–Mouse Orthologue Pairs. We identified 9,207 genes whose human and mouse official gene symbols were identical. We then selected only the nuclear protein-coding genes, which amounted to 7,408 genes. We regarded these gene pairs as orthologues. Among these, we were able to collect 9-kb genomic upstream nucleotide sequences without any described genic regions for 3,750 orthologous gene pairs. Then, we excised 8-kb stretches upstream of the translation start sites. We set a 1-kb margin to decrease the frequency of cases in which the excised 8-kb sequences overlapped with promoter regions or 3′ regulatory sequences of adjacent genes. For all of the 3,750 pairs of human and mouse genes, the accessions of the contig entries used are shown in Table 3, which is published as supporting information on the PNAS web site, together with the positions of the excised sequences. Finally, we were able to make a global pairwise alignment for each of the 3,750 orthologue pairs by using realigner.

Simulation Analysis. To estimate the level of the number of identical sites obtained by chance in the 8-kb alignments, we conducted a simulation study in which 10,000 pairs of randomly generated 8-kb sequences were processed in the same way as described above (Fig. 1). Because 10 of the 10,000 random sequence pairs exceeded 287 as the number of identical sites, we regarded any orthologous genomic global alignment whose number of identical sites was below this value as not significant at the P < 0.001 level and only used the orthologues above this threshold for further analyses (3,055 of 3,750 orthologue pairs).

Validation for Conserved Sequences Reflecting Regulatory Regions. To confirm that the conserved upstream sequences obtained by the alignment procedures were actually able to reflect the major regulatory regions, we focused on the experimentally validated promoter regions that were registered in the EPD. We could assign 266 entries in the EPD to the human upstream sequences under analysis. As a result, 80.0% of the 266 EPD sequences were covered by conserved sequences in the alignments over at least a 10-bp stretch, and, for 40.0% of the 266 entries, the entire 60-bp EPD promoter region was completely covered in the conserved sequences. These results indicate that our genomic alignment procedures are effective for detecting the regulatory regions.

Repetitive Sequences of Highly Upstream-Conserved Genes. Within all of the upstream sequences of the top 1–10 upstream-conserved orthologue pairs, we found very few transposable element (TE) insertions (i.e., 0.42% of the human and 0.56% of the mouse 8-kb upstream sequences) by using repeatmasker (http://ftp.genome.washington.edu/RM/RepeatMasker.html), although TE insertions have been reported to occupy ≈45% of the human (13) and 37.5% of the mouse (7) genome sequences. This result indicates that the highly conserved upstream sequences that we identified have been under strong selection against TE insertions.

High Degree of Upstream Conservation for TF Genes. We have listed the top 30 upstream-conserved orthologues in Table 1 (see also Fig. 2 and the list of all of the 3,055 genes under analysis shown in Table 4, which is published as supporting information on the PNAS web site). Notably, we found that nine of the 10 most upstream-conserved orthologues were TF genes (P < 2 × 10^–8). Within the top 30 upstream-conserved genes, 19 were also TF genes (P < 5 × 10^–12) (Table 1). The high occurrence of TF genes also continued for the top 200 upstream-conserved genes, i.e., 62 of the top 200 genes encoded TFs (P < 5 × 10^–15), although TF genes occupy only 10.0% (305/3,055) of the genes analyzed.

Relationship Between the Degree of Upstream Sequence Conservation and Gene Function. To further elucidate the relationship between the gene function and the degree of upstream sequence conservation, we subcategorized the genes by their function, based on the GO annotation, into the following six categories: (i) metabolism, (ii) biosynthesis, (iii) cell cycle, (iv) signal transduction, (v) developmental process, and (vi) TFs. As a result, we found that the TF genes and developmental process-related genes showed the highest degrees of upstream sequence conservation within the six categories (Fig. 3). In more detail, 61 developmental process-related genes appeared within the top-200 upstream-conserved genes (P < 5 × 10^–12), although such genes occupy only 11.4% (347/3,055) of the genes in the present analysis.

In contrast to the TF genes and developmental process-related genes, genes involved in metabolism and biosynthesis showed less upstream sequence conservation. These categories of genes belonged more frequently to the lower ranks of upstream sequence conservation. An intermediate degree of upstream conservation was shown for cell cycle- and signal transduction-related genes, which were distributed almost evenly throughout all of the ranks concerning the degree of upstream conservation (Fig. 3). Although a considerable number of cell cycle-related genes contribute to common processes among eukaryotes, the degree of upstream sequence conservation remained much less than that of the developmental process-related genes and TF genes. Our results clearly show that the TF genes have a remarkably high degree of conservation in the upstream regions that are considered to play major regulatory roles.

Higher Conservation of Upstream Regions of Developmental Process-Related TF Genes. We further focused on the genes that were involved in both TF activity and developmental processes. These developmental process-related TF genes gave an even higher degree of upstream sequence conservation. Notably, >27% of the developmental process-related TF genes appeared in the top 200 upstream-conserved rank. They also appeared more frequently in the top three ranks (i.e., the top 1–800 upstream-conserved genes) than the TF genes that are not involved in developmental processes (P < 1 × 10^–6; rank-sum test) and non-TF genes that are involved in developmental processes (P < 1 × 10^–15; rank-sum test) (Fig. 4). We showed that developmental process-related TF genes have a much higher degree of upstream sequence conservation.

Candidates for Regulatory SNPs in the 8-kb Upstream-Conserved Sequence of the ZFHX1B Gene. According to the RefSeq contig annotation, we found that 11 SNPs were located within the 8-kb stretch upstream of the human ZFHX1B gene, which was the most upstream-conserved orthologue between humans and mice. Of these, six were within the human-mouse conserved sequences that we identified (Table 2). We showed the physical locations of all of the conserved sequences (96,873 sequence blocks) that we identified in the alignments of the 3,055 human–mouse orthologous 8-kb upstream sequence pairs in Table 5, which is published as supporting information on the PNAS web site.

Discussion

High Upstream Conservation for Development-Related TF Genes. In the present study, we have shown that the genes with high upstream sequence conservation are predominantly TF genes. More interestingly, we have also elucidated that developmental process-related TF genes have an even higher degree of upstream sequence conservation.

Because our resultant upstream conserved sequences could reflect regulatory sequences, these TF genes may have a larger number of regulatory sequences and may form highly connected regulatory networks between TF genes, such as those found in yeast by using genome-wide experiments (14). Because such highly connected biological networks have been reported to be preferentially retained evolutionarily (15), the potentially high connectivity of the regulatory networks that are constituted of TF genes (in particular, developmental process-related TF genes) might contribute to the evolutionary stasis.

Therefore, the surplus degree of upstream conservation we found for the developmental process-related TF genes can illustrate the evolutionary constraint that has been considered to be inherent in the program of ontogeny, or so-called developmental constraint (16, 17).

Pleiotropism of ZFHX1B Could Explain the Requirement of Its Strict Regulation. We found that the most upstream-conserved human–mouse orthologous gene, ZFHX1B, had 26 conserved sequence blocks and showed 75% nucleotide identity on average along the 8-kb upstream sequence. The ZFHX1B gene encodes a TF, Smad-interacting protein 1 (SIP1) (18), and mutations in the exons of this gene have been reported to cause a form of Hirschsprung disease (19, 20), the most common congenital malformation (≈1 in 5,000 live births). Besides Hirschsprung disease, mutations in one allele of the ZFHX1B gene have been shown to cause a wide spectrum of congenital anomalies (19, 20), such as microcephaly and malformations of the brain (cerebral atrophy or callosal agenesis), heart (patent ductus arteriosus), face morphology, and genitourinary organs (21), along with complex developmental disorders, such as mental retardation, delayed motor development, and epilepsy (19–21). These findings indicate that haploinsufficiency of the ZFHX1B gene can have a variety of influences on multiple organs in multiple phases. This wide breadth of pleiotropism of the ZFHX1B gene (20) should require strict temporal and spatial control of its expression, which could account for the large number of upstream regulatory sequences.

Candidates for Regulatory SNPs for Further Examination. Although all of the mutations in ZFHX1B reported to cause disease phenotypes have been within the exons, it would also be meaningful to examine nucleotide changes in the highly conserved noncoding upstream sequences of ZFHX1B that we identified, because partial dysregulation of this gene could cause milder forms of those clinical features. Alternatively, such mutations could be genetic factors relevant to other complex diseases, because of its wide breadth of pleiotropism.

We found 11 SNPs located within the 8-kb stretch upstream of ZFHX1B, six of which were included in the conserved sequences (Table 2). These six SNPs could be candidates for regulatory SNPs. Because we identified 96,873 conserved sequence blocks in the 3,055 pairs of human–mouse orthologous 8-kb upstream sequences (shown in Table 5), these evolutionarily conserved sequences may be good targets for identifying novel mutations that are related to diseases or disease susceptibilities, and they may also be useful as targets for identifying regulatory sequences in high-throughput experimental analyses.

In summary, the present study has demonstrated the possible molecular basis of evolutionary constraint related to the regulatory program of ontogeny as the surplus degree of upstream sequence conservation for developmental process-related TF genes. By means of cross-species comparison of genome sequences and the utilization of functional annotations, we have demonstrated the gene expression regulatory features characteristic of TF genes and development-related genes, which await further experimental evidence on a genome-wide scale.