Whole-exome sequencing analysis in 10 families of sporadic microtia with thoracic deformities.

BACKGROUND: Microtia is a congenital malformation of the external ear and may occur as an isolated deformity or as part of a syndrome. Our previous study found a high correlation between microtia and thoracic deformities, thus, we propose that external ear and thorax development may be regulated by certain genes in common.
METHODS: We performed exome sequencing on 10 families of sporadic microtia with thoracic abnormalities. We identified mutated genes under different models of inheritance, and checked them through Mouse Genome Informatics and association analysis.
RESULTS: We identified 45 rare mutations, including 9 de novo mutations, 20 heterozygous mutations, 3 homozygous mutations, and 13 hemizygous mutations, of which 2 are likely to be causative. They are de novo missense variant in PHF5A and compound heterozygous mutations in CYP26B1, of which CYP26B1 mutation is highly likely pathogenic.
CONCLUSION: The results indicate that certain genes may affect both external ear and thorax development, and demonstrate the benefits of whole-exome sequencing in identifying candidate genes of microtia. This study provides a new way for genetic exploration in microtia.

INTRODUCTION

Microtia is a congenital anomaly of the external ear that ranges in severity from mild structural deformities to complete absence of the auricle with stenosis or atresia of the external auditory canal (Suutarla et al., 2007). The prevalence rates reported vary from 0.83 to 4.34 per 10,000 births (Cox et al., 2014). Microtia may occur as an isolated deformity or as part of a spectrum of anomalies or a syndrome. Roughly 20 to 60% of microtia patients have associated deformities or an identifiable syndrome (Hartzell & Chinnadurai, 2018; Shibazaki‐Yorozuya & Nagata, 2019; Stoll et al., 2016). The common deformities include: cleft lip and/or palate, micrognathia, mandibular dysplasia, macrostomia, facial asymmetry, renal abnormalities, cardiac defects, and polydactyly. The severity of microtia was classified into grades I, II, or III according to the classification used by Weerda, (1988). In brief, grade I microtia exhibits only mild deformity, with the auricle being slightly smaller than normal, each part of which can be clearly distinguished. In grade II microtia, the size of the auricle is one‐half to two‐thirds of the normal size and its structure is only partially retained. In grade III microtia, the auricle is severely malformed and usually exhibits a peanut shape. The extreme case where there is no external ear and auditory canal is called anotia or microtia grade IV.

In our previous study, we found the incidence of thoracic deformities, including rib, rib cartilage, and spine deformities, was high in patients with microtia. And the incidence of thoracic deformities was correlated with the grade of microtia: the poorer one auricle developed, the higher the incidence of thoracic deformities (Wu et al., 2015; Yang et al., 2016). It has been reported that rib deformities may present in some syndromes, for example, thoracic outlet syndrome (Jones et al., 2019), Jeune syndrome (asphyxiating thoracic dystrophy, ATD; Page et al., 2017), and Poland syndrome (Romanini et al., 2016). Rib cartilage anomalies have been rarely reported. Spinal deformities may occur in isolation, as well as in association with multisystem malformations such as oculo‐auriculo‐vertebral spectrum (Beleza‐Meireles et al., 2015; Renkema et al., 2019; OAVS). In order to better understand the association between microtia and thoracic deformities, the study just focuses on those with microtia and thoracic involvement and no other feature. The embryological development of the ribs and vertebrae is closely associated, and the complete anatomical pattern is formed in the mesenchyme during the first 6 weeks of intrauterine life (Tsou et al., 1980). Axial skeleton and cartilage formation in the embryo both involve the progeny of mesenchymal stem cells. Therefore, microtia and thoracic anomalies are likely to be the result of a common pathogenetic mechanism occurring during embryonic development. We postulate that the thoracic deformities may represent a new syndrome characterized by “microtia‐thorax anomaly” or another extended manifestation of microtia, which suggests that malformations of the external ear and thorax may reflect the actions of certain genes in common. Further study examining or selecting responsible genes in eligible microtia cases with thoracic deformities are needed to investigate the pathogenetic mechanism involved.

To date, few studies have focused on microtia patients with thoracic deformities and their candidate genes, although some studies have identified causal genes for microtia and related syndromes. In recent years, the technology of whole‐exome sequencing (WES) has proven to be effective and feasible and has been successfully used in identifying the cause of many diseases (Iglesias et al., 2014; Linthorst & Hollak, 2019; Sun et al., 2019), which also offers a new platform for exploring the candidate genes of microtia and have made some progress. Brown et al., (2013), sequenced the whole exome of a three‐generation microtia family and defined HOXA2 (OMIM *604685) haploinsufficiency as the first genetic cause for autosomal‐dominant nonsyndromic microtia. By WES, Lopez et al. (2016), detected a heterozygous mutation in the myelin transcription factor 1 (MYT1, OMIM *600379) gene in one OAVS patient, and functional studies identified MYT1 as the first gene implicated in OAVS. Vetro et al., (2017), performed exome sequencing in a patient with a clinical diagnosis of Meier–Gorlin syndrome (characterized by microtia, primordial dwarfism, and patellar aplasia/hypoplasia) and identified MCM5 (OMIM *602696) as a novel causative gene. Here, we performed exome sequencing on a cohort of 10 sporadic microtia patients (along with their parents) with diverse thoracic abnormalities. Through bioinformatic analysis of the sequencing data, we attempted to identify genes or mutations that may be associated with microtia and thoracic deformities.

MATERIALS AND METHODS

Ethical compliance

This study was approved by the institutional ethics committee of Plastic Surgery Hospital, Chinese Academy of Medical Sciences and Peking Union Medical College, and written informed consent was obtained from all subjects and their parents.

Recruitment of subjects

Ten microtia patients (four females, six males; mean age, 7.8 years; range, 5–12 years) and their nonconsanguineous parents were recruited from Plastic Surgery Hospital, Chinese Academy of Medical Sciences and Peking Union Medical College. All patients were asked to provide a detailed family history, and none had a family history of external ear malformation. Ten patients all had thoracic deformities, including eight cases with 12th‐rib anomalies, one case with first and second ribs anomalies, and one case with 12th‐rib as well as rib cartilage anomalies (Table 1). Patients with chromosome abnormality or syndromic microtia such as Treacher Collins syndrome, OAVS, CHARGE syndrome, branchio‐oto renal (BOR) syndrome, or Jeune syndrome, were excluded from the study.

TABLE 1

Clinical characteristics of 10 patients

Sample ID	Sex	Age (years)	Side of microtia	Thoracic anomalies
S1	Female	7	Right	Bilateral absent ribs of 12th
S2	Female	12	Left	Bilateral absent ribs of 12th
S3	Male	10	Right	Right hypoplastic rib of 12th
S4	Male	5	Left	Left hypoplastic rib of 12th
S5	Female	7	Right	Bilateral absent ribs of 12th
S6	Male	5	Right	Bilateral hypoplastic rib of 12th
S7	Male	8	Left	Bilateral hypoplastic rib of 12th
S8	Female	11	Left	Bilateral absent ribs of 12th
S9	Male	7	Right	Right hypoplastic ribs of first and second
S10	Male	6	Right	Bilateral absent ribs of 12th, fused rib cartilages of left second and third

Whole‐exome sequencing

For detailed methodology see Appendix S1. In brief, patients’ and their parents’ DNA samples were exome sequenced using Agilent liquid capture system (Agilent SureSelect Human All Exon V5) according to the manufacturer's protocol followed by paired‐end sequencing (150 bp reads) on Illumina HiSeq 4000. After quality control, high‐quality clean data were obtained and mapped to the reference genome (UCSC hg19) by Burrows–Wheeler Aligner (BWA) software to get the original mapping result in BAM format. Subsequently, Samtools and Picard (http://broadinstitute.github.io/picard) were spectively utilized to sort bam files, do duplicate‐marking to generate final bam file. Samtools mpileup and bcftools were used to do variant calling and identify SNP, indels. ANNOVAR and a variety of databases were applied for functional annotation of variants.

Identification of candidate genes or variants

First, to identify the most likely pathogenic de novo mutations, variants obtained from previous steps were then filtered with the minor allele frequency (MAF) ≥0.5% in 1000 Genomes Project Consortium (1000 Genomes). Only SNVs occurring in exons or in canonical splice sites (splicing junction 10 bp) were further analyzed. Synonymous SNVs which did not lead to an amino acid alternation were discarded to get nonsynonymous SNVs. Then, the retained nonsynonymous SNVs were submitted to PolyPhen, SIFT, Mutation Taster, Combined Annotation‐Dependent Depletion (CADD) for functional prediction. The SNVs that were predicted to be not benign by at least two softwares could be retained. To identify inherited recessive and X‐linked SNVs and indels, we first merged the variant call format (VCF) files from the three variant callers for each individual of the family, and then, identified genes harboring rare, high quality, functional variants (predicted protein consequences were essential splice site, stop gained, frameshift coding, non‐synonymous, and stop lost) under recessive and X‐linked models.

Second, to further determine the mutated genes involved in microtia, we checked the rare mutated genes identified under different genetic type with Mouse Genome Informatics (MGI, http://www.informatics.jax.org/). Genes that expressed in branchial arches or influenced the proliferation, migration, and differentiation of cranial neural crest cells (CNCC) were retained.

Last, to test whether the above genes were related with clinical genotypes, we conducted a phenotype‐based prioritization of candidate genes involved in microtia and rib abnormalities using a computational tool called phenolyzer (http://phenolyzer.usc.edu). By entering microtia and rib abnormalities respectively as search name, prioritized gene list was obtained, and then, applied to screen out the candidate genes and track their information on gene–gene interaction.

RESULTS

Whole‐exome sequencing in 10 sporadic microtia patients with thoracic deformities, along with their parents, was performed (a total of 30 individuals). The average sequencing depth on targeted coding regions was 123×. A mean of 99% of bases in the targeted coding regions were covered by at least 10 reads. Number of bases with a phred‐like calibrated quality score of 20 or above was greater than 90%. Number of bases with quality score of 30 or above was greater than 80% (Figure S1 and Table S1). Mean error rate was lower than 0.1%, demonstrating the high quality of the sequencing. No family history of auricular abnormalities were provided, suggesting dominant inheritance is unlikely. Therefore, we evaluated different classes of potentially pathogenic, rare coding variants under dominant de novo, recessive and X‐linked modes of inheritance.

Through a strict filter, we identified 45 rare coding mutations, including 9 de novo mutations, 20 heterozygous mutations, 3 homozygous mutations, and 13 hemizygous mutations from 10 microtia patients (Table 2). Nine de novo SNVs and indels were classified as follows: six missense mutations, one frameshift insertion, one stop gain mutation, and one unknown mutation (Table S2). The six missense mutations included c.G418A (p.D140 N) in HDAC11 (OMIM *607226, NM_001136041.3), c.C12151G (p.P4051A) in CMYA5 (OMIM *612193, NM_153610.5), c.G2522A (p.R841H) in CTC1 (OMIM *613129, NM_025099.6), c.G3892A (p.G1298R) in GON4L (OMIM *610393, NM_001282856.1), c.A424G (p.I142V) in KIN (OMIM *601702, NM_012311.4), and c.G904 T (p.D302Y) in SMG6 (OMIM *610963, NM_017575.5). One frameshift insertion was c.326_327insGTGTCTCG (p.R109fs) in VWA7 (OMIM *609693, NM_025258.3). One mutation, c.C162A (p.Y54*) in PHF5A (OMIM *617846, NM_032758.4), was predicted to lead to a truncated protein owing to a premature stop codon. One unknown mutation was found in MUC19 (OMIM *612170, NM_173600.2).

TABLE 2

Candidate genes identified in 10 microtia patients with thoracic deformities

Sample ID	Genes with de novo mutations	Genes with inherited autosomal mutations (recessive, compound heterozygous)	Genes with inherited autosomal mutations (recessive, homozygous)	Genes with inherited mutations on X chromosome hemizygous)
S1	HDAC11	ADGRV1; TBXAS1	OR4X1, ZNF2	—
S2	PHF5A	OGDH; VLDLR	AKAP12	—
S3	—	CELSR1	—	PHKA1; PLXNB3
S4	VWA7	CYP26B1; RSPH6A; TTYH1	—	GAB3; PLXNA3; RGAG1; RRAGB
S5	MUC19	SLC4A4; XIRP2; ZNF665	—	—
S6	—	ERICH3; KIAA2012		CDKL5; STARD8
S7	CMYA5; CTC1; GON4L	AHNAK; KIAA1217; PNPLA7; TTN	—	BEND2; L1CAM; ZNF75D
S8	KIN; SMG6	MUC16	—	—
S9	—	—	—	PDHA1
S10	—	CHST13; PTPRS	—	PFKFB1

In recessive and X‐linked inheritance pattern, variants (SNVs and indels) were identified by filtering for rare (MAF < 0.5%) and functional variants. We identified a mean of 3.6 mutated genes per patient (range of 0–4) with a total of 36 different mutated genes across the 10 patients, containing 51 rare functional variation sites and three unannotated variation sites. No common mutations were found in any two patients. Of all the functional variants, 46 were missense, 3 were stopgain, 1 was frameshift deletion, and 1 was nonframeshift deletion (Table S3).

Next, all the mutated genes were checked by MGI. For the de novo mutations, we found two genes, PHF5A and KIN, expressed in branchial arches (Figures S2–S3). And based on existing phenotypes in knockout mice, no genes were found to have an influence on the proliferation, migration, and differentiation of CNCC. For the inherited recessive and X‐linked mutations, AKAP12 (OMIM *604698, NM_005100.4), CYP26B1 (OMIM *605207, NM_001277742.1), and AHNAK (OMIM *103390, NM_001620.3) were detected to express in branchial arches, and no genes were found to cause the phenotypes related to the development of CNCC (Figures S4–S6).

At last, we utilized phenolyzer to detect whether the genes identified through the above steps were related with clinical genotypes. When microtia, as a query term, was input in the interface of phenolyzer, we got 10348 candidate genes, including PHF5A ranked 6627 and CYP26B1 ranked 1516 (Tables S4–S5). Interestingly, when rib abnormalities were input in the interface of phenolyzer, a total of 12242 candidate genes were provided, of which PHF5A ranked 8073 and CYP26B1 ranked 1012 (Tables S6–S7).

The de novo variant in PHF5A was found in S2, a girl with left microtia and bilateral absent 12th ribs. It resulted in substitution of a alanine for a cytosine at position 162 (C162A), which led to a premature stop codon. The variant was neither reported in 1000 Genomes, Exome Aggregation Consortium (ExAC), ESP6500 databases nor recorded in Genome‐Wide Association Study (GWAS) and Human Gene Mutation Database (HGMD). The variant was also not found in Genome Aggregation Database (gnomAD), and the pLI score in gnomAD is 8.6605e‐01. This stop gain mutation was predicted to be not only “disease_causing_automatic” by Mutation Taster, but also conservative by SiPhy (score 10.8538), PhyloP (score 3.786), and GERP (score 4.47). Besides, the Phred‐scaled score of CADD was 38, suggesting the mutation was likely to be pathogenic.

In S4, a male microtia patient with left hypoplastic rib of 12th, the compound heterozygous variations detected in CYP26B1, c.G476A (p.R159Q), and c.G370A (p.G124S), were inherited from the mother and the father, respectively. The allele frequencies of both variants were very low in ExAC, 0.0000417 and 0.00000825, respectively. The two variants neither appeared in 1000 Genomes and ESP6500 databases, nor reported in GWAS and HGMD. The gnomAD allele frequencies were also very low, 0.00004719 and 0.00000812, respectively. Almost all of the occurrences have been in East Asians, which is consistent with our study subjects. The result may suggest the variants of CYP26B1 are more meaningful to East Asians. In silico pathogenicity prediction was conducted on the two missense variants and the results from at least two softwares indicated them as “deleterious/probably damaging/disease_causing.” Above all, according to the phenolyzer results, CYP26B1 is closely related to phenotypes and hence highly likely pathogenic. Subsequently, we analyzed the sequences of CYP26B1 and its coding protein, finding that the two variation sites were both in conserved domain. And protein coding sequences from higher mammals were aligned to demonstrate that the two sites were both highly conserved (Figures 1 and 2). And multiple alignments of 100 vertebrate species using the BLAT search on the UCSC Genome Browser (https://genome.ucsc.edu/) also showed the two sites were highly conserved. The conformational space of CYP26B1 coding protein was predicted by SWISS‐MODEL (https://www.swissmodel.expasy.org/). The three‐dimensional structures showed no difference between mutated protein and normal protein, except for some differences in side‐chain groups.

FIGURE 1

Protein coding sequences from different species demonstrated that the site 159 in CYP26B1 was highly conserved

FIGURE 2

Protein coding sequences from different species demonstrated that the site 124 in CYP26B1 was highly conserved

DISCUSSION

Microtia is a common deformity of embryonic development seen in the plastic surgery clinic which may present along with other malformations, therefore, individuals with microtia should be examined for other dysmorphic features. In our preliminary clinical study (Yang et al., 2016), a total of 214 microtia patients received a preoperative three‐dimensional chest computed tomography and we found that the incidence of thoracic deformities (including rib, rib cartilage, and spine deformities) was higher in microtia patients and positively correlated with the grade of microtia. That is, from the aspect of incidence, patients with microtia III were observed to have a higher incidence of thoracic deformities than those with microtia II. Patients with microtia II were found to have a higher incidence of rib and spinal deformities than those with microtia I. The high incidence of thoracic deformities and their co‐occurrence with microtia suggests that there may have been a common pathogenetic mechanism at the same stage of embryologic development. Thus, we propose that external ear and thorax development may be regulated by certain genes in common. To date, no studies have reported an investigation of the associations between external ear malformations and thoracic deformities etiology.

Microtia is a clinically heterogeneous spectrum and thus many more genes and causal variants are likely to be identified. A variety of genetic methods were used to study microtia, the four most common including GWAS, WES, linkage studies in large families, and copy number variation investigations. Especially, the WES has been considered to improve the discovery of new candidate genes and the sensitivity of detecting de novo and compound heterozygous variants (Bekheirnia et al., 2017; Qiao et al., 2019; Sohn et al., 2018; Yang, Shen, et al., 2019).

In this study, we first performed exome sequencing of 10 sporadic microtia patients (along with their parents) with diverse thoracic deformities and identified nine de novo functional mutations and 36 inherited functional mutations, including 20 heterozygous mutations, 3 homozygous mutations, and 13 hemizygous mutations. Second, all the mutated genes were screened by the MGI database. Five genes, PHF5A, KIN, AKAP12, CYP26B1, and AHNAK, were discovered to be expressed in branchial arches. Last, we used phenolyzer to detect the association between phenotypes and above genes, and then, screened out two candidate genes, PHF5A and CYP26B1, which are considered to be causal.

PHF5A belongs to a small transcription factor or co‐factor which encodes 110 amino acids. It is found in many different species ranging from yeast to human, playing an important role in RNA splicing, cell cycle, and DNA damage repair processes. Deletion of PHF5A in yeast is lethal and knockdown in C. elegans results in aberrant organogenesis during early development, suggesting its importance for embryo formation and tissue morphogenesis (Trappe, Ahmed, et al., 2002; Trappe, Schulze, et al., 2002). By reviewing the related literatures, we have found most studies focus on the role of PHF5A in cancer so far, and have showed that knockdown of PHF5A significantly suppresses cell proliferation, migration, and invasion of tumor (Yang, Zhang, et al., 2019; Zheng et al., 2018). Another study demonstrated that PHF5A was essential for maintaining pluripotency, cellular reprogramming, and myoblast specification (Strikoudis et al., 2016).

According to the results of phenolyzer, PHF5A is related to not only microtia, but also rib abnormalities. It is worth mentioning that we find the same biosystem in the association analysis of PHF5A and “microtia” as well as “rib abnormalities.” That is, PHF5A is in the same gene expression with the seed gene, SMAD4 (OMIM *600993, Table S2 and S4). SMAD4 is an essential mediator of transforming growth factor‐β/bone morphogenetic protein (TGF‐β/BMP) signaling, widely expressed in many tissues, for example, cartilage growth plate, skin, and heart. To the best of our knowledge, the specific role of SMAD4 during ear development has been studied in two studies so far, which demonstrate that SMAD4 is essential for proper formation of outer ear cartilage, inner ear development, and normal auditory function (Teng et al., 2011; Yang et al., 2009). SMAD4 mutant mice exhibit a short ear phenotype due to defective chondrocyte maturation and cartilage production (Teng et al., 2011), which suggests that SMAD4 is important for cartilage development of the outer ear. Thus, from the above, we assume that de novo missense variant in PHF5A is possibly pathogenic. However, there have been no reports or animal models about PHF5A mutations leading to microtia or thoracic deformities, which requires further confirmatory genetic and functional studies.

Additionally, the compound heterozygous mutations in CYP26B1 is proposed to be highly likely pathogenic. CYP26B1 is one of the three CYP26 gene isoforms (CYP26A1, CYP26B1, and CYP26C1) which encode the cytochrome‐P450 enzymes that catabolize retinoic acid (RA). RA is the active metabolite of vitamin A, acting as an essential molecule of cell–cell signaling during development of vertebrate embryos (Duester, 2008). Too little and too much RA both cause the associated malformation syndromes, involving ears, eyes, central nervous, muscle, skeleton, and urogenital system (Lee et al., 2012; Petrelli et al., 2019). In the zebrafish, CYP26B1 deletion has been shown to result in overall defective craniofacial cartilage development with smaller head, severely decreased number of vagal branchiomotor neurons and defective or absent jaw cartilage (Reijntjes et al., 2007). Similarly, CYP26B1 null mice exhibit craniofacial abnormalities, including small external pinnae, cleft palate, micrognathia, limb and other bone, and cartilage deformities (Maclean et al., 2009; Minegishi et al., 2014). In addition, deletions of CYP26B1 gene have also been reported in humans. The previous work of Laue et al. (Laue et al., 2011), pointed toward homozygous mutations of CYP26B1 as being responsible for the generation of lethality, skeletal, and craniofacial abnormalities. The surviving individuals carrying a CYP26B1 deficiency may present a large number of developmental abnormalities, such as intellectual disability, facial asymmetry, ear malformations, skeletal dysplasia, vertebral abnormalities (Li et al., 2018; Laue et al., 2011; Morton et al., 2016). And further study showed that mutations in CYP26B1 reduced the catalytic activity of the encoded enzyme to near zero or 30%, implying that localized accumulation of RA could have deleterious effects on morphogenesis, most especially in the developing skeleton. Hence, based on these reports and animal models, we conclude that CYP26B1 mutation is highly likely pathogenic, although p.R159Q and p.G124S have not previously been reported. Furthermore, it is worth mentioning that CYP26B1, as the same as PHF5A, is also in the same gene expression with SMAD4. This finding may indicate the overlapping pathway/network of different genes, but further research is needed.

In short, our study identifies two candidate genes, PHF5A and CYP26B1, which may be associated with microtia and thoracic deformities. Especially, CYP26B1 is the most plausible gene. Therefore, our next step is to verify its function and establish an animal model if necessary. The evidence in favor of PHF5A is relatively weak, so we still need to detect mutation of PHF5A in sporadic microtia patients in the following work. If mutation of PHF5A is detected in another microtia patient, gene function experiment will be performed. In conclusion, this is the first study to demonstrate the genetic landscape in sporadic microtia patients with thoracic deformities. Whole‐exome sequencing is an effective strategy to elucidate the etiology of congenital diseases, and our findings provide a basis for further genetic research on microtia with thoracic deformities and broaden our mind to explore the pathological mechanisms of microtia.

CONFLICT OF INTEREST

The authors declare no conflicts of interest.

AUTHOR CONTRIBUTIONS

HY and BP designed the study. MY, YZ, and CW did the clinical assessment and recruited the patients and their parents. MY wrote the original draft. All authors contributed to data analysis and revising the manuscript. All the authors gave approval of the final version.

Supplementary Material

Click here for additional data file.