Mandibulofacial Dysostosis Attributed to a Recessive Mutation of CYP26C1 in Hereford Cattle

In spring 2020, six Hereford calves presented with congenital facial deformities attributed to a condition we termed mandibulofacial dysostosis (MD). Affected calves shared hallmark features of a variably shortened and/or asymmetric lower mandible and bilateral skin tags present 2–10 cm caudal to the commissure of the lips. Pedigree analysis revealed a single common ancestor shared by the sire and dam of each affected calf. Whole-genome sequencing (WGS) of 20 animals led to the discovery of a variant (Chr26 g. 14404993T>C) in Exon 3 of CYP26C1 associated with MD. This missense mutation (p.L188P), is located in an α helix of the protein, which the identified amino acid substitution is predicted to break. The implication of this mutation was further validated through genotyping 2 additional affected calves, 760 other Herefords, and by evaluation of available WGS data from over 2500 other individuals. Only the affected individuals were homozygous for the variant and all heterozygotes had at least one pedigree tie to the suspect founder. CYP26C1 plays a vital role in tissue-specific regulation of retinoic acid (RA) during embryonic development. Dysregulation of RA can result in teratogenesis by altering the endothelin-1 signaling pathway affecting the expression of Dlx genes, critical to mandibulofacial development. We postulate that this recessive missense mutation in CYP26C1 impacts the catalytic activity of the encoded enzyme, leading to excess RA resulting in the observed MD phenotype.


Introduction
Over 250 Mendelian traits in cattle are reported in the Online Mendelian Inheritance in Animals database (https://omia.org/home). Often deleterious syndromes in cattle are attributed to variants inherited in an autosomal recessive manner [1]. Due to this inheritance pattern, clinical signs of disease may not appear for many generations after the causal mutation originates. However, artificial selection in livestock and the commonplace use of artificial insemination and embryo transfer can expedite widespread proliferation of a deleterious variant. Once a deleterious defect is identified, prompt identification of carrier animals is necessary to prevent economic loss. The significant impact of a single deleterious variant in livestock can be exemplified by a recessively inherited mutation in APAF1, traced back to a Holstein bull born in 1962. The recessive genotype, detrimental to cow fertility, was estimated to have resulted in 525,000 abortions costing the industry $420 million [2].
In March and April of 2020, three herds reported a total of six purebred Hereford calves born with unusual defects of the face and jaw attributed to a condition we termed mandibulofacial dysostosis (MD). One of the reporting herds also noted a calf born with a similar phenotype in 2018. The calves were live born, of normal size, and with an intact suckle reflex but weak suckling ability. Grossly, several calves appeared to have a widened upturned "smile," a variably shortened and/or asymmetric lower mandible, and unique, bilateral skin tags just caudal to the commissure of the lips. The three herds reporting MD cases were each in a different state (Iowa, Wyoming, Missouri) with typical summer grazing and winter feeding programs, making an environmental cause unlikely. Sires and dams of affected calves were consanguineous.
The similarity in the reported clinical description among the affected calves, pedigree analysis, and description of similar phenotypes in children [3] and mice [4,5] led to the hypothesis that a de novo, autosomal recessive mutation may be causative of this novel condition in Hereford cattle. Pathologic findings of retained Meckel's cartilage in the affected calves further suggested such a mutation might disrupt development of the first pharyngeal arch (PA1), possibly through endothelin-1 (ET1) and Dlx signaling. Expression of Dlx homeobox genes in the cranial neural crest cells (CNCC) of the embryonic PA1 provides patterning information during jaw formation [4,6]. More specifically, ET1 signaling in PA1 activates Dlx5 and Dlx6, driving differentiation of the mandible and maxilla [5,6]. The identification of a causative variant would allow breeders to identify carrier individuals to avoid production of affected calves; it also has the potential to yield novel information regarding the regulation of craniofacial development. To investigate our hypothesis, multiple affected calves were obtained to establish a phenotypic characterization of the defect and case definition. DNA from these calves, family members, and herd mates was collected for whole-genome sequencing (WGS) to identify possible de novo variation impacting mandibulofacial development.

Materials and Methods
Animals utilized in this study were sampled in compliance with approved University of Nebraska IACUC project number 1970: Diagnostic Investigation into Natural Animal Disease Events. Five affected calves, two heifers and three bulls, were received at the University of Nebraska-Lincoln Veterinary Diagnostic Center (UNL-VDC) for autopsy evaluation, four after euthanasia on the farm and one live, within 24 h of birth. The live calf was evaluated for hearing and vision and then euthanized by intravenous injection of pentobarbital sodium and phenytoin sodium (Euthasol, Virbac AH, Inc., Fort Worth, TX, USA).
Sire and dam of all reported, affected calves as well as their extended pedigrees (records available through the American Hereford Association, https://hereford.org) were evaluated to identify common ancestors. The herd of origin and date of birth was also noted.
Tissue samples (ear) were collected from the 5 affected calves received at the UNL-VDC. Semen or whole blood samples were obtained from parents of affected calves (as available) and initially from 7 other related individuals. All samples were stored at −20 • C. DNA was isolated from EDTA blood samples with a Gentra Puregene Blood Kit (Qiagen, Venlo, Netherlands) utilizing the following modified protocol. To obtain the buffy coat, blood tubes were centrifuged at 2000× g (15 min, 4 • C). 250 µL of buffy coat was combined with 900 µL red blood cell lysis solution, vortexed, and incubated (5 min, 22 • C). Samples were then centrifuged (13,000× g, 2 min) and the supernatant discarded. 450 µL of red blood cell lysis solution was added to the pellet, vortexed, and incubated (5 min, 22 • C). Samples were again centrifuged (13,000× g, 2 min) and supernatant discarded before adding 900 µL cell lysis solution and 6 µL Proteinase K. The samples were vortexed then incubated (30 min, 35 • C). After cooling on ice to room temperature, 200 µL of protein precipitation solution was added, the samples vortexed, and incubated on ice (5 min). The tubes were centrifuged (13,000× g, 2 min). To precipitate DNA, the supernatant was poured into a new tube containing 800 µL of 100% isopropanol, the tubes were inverted 50 times, and centrifuged (8000× g, 2 min). The supernatant was discarded and pellet dried for 1 min before washing with 300 µL of 70 % ethanol. The cleaned pellet was dried for 15 min before the DNA was rehydrated in 100 µL of DNA hydration solution overnight (22 • C) before storage at 4 • C.
Filtering to identify candidate variants was completed utilizing VCFtools (https://vcftools.github. io). Candidate variants were identified as those homozygous for the alternative allele in all affected calves, found in a heterozygous state in all obligate carriers (parents of affected calves), and either heterozygous or homozygous for the reference allele in the related individuals. Variants were further prioritized by evaluating their predicted impact on gene/protein function as annotated by the Ensembl Variant Effect Predictor (https://uswest.ensembl.org/info/docs/tools/vep/index.html) and SNPEff [11].
Sequence Read Archive (SRA) data were obtained for animals mapped to the ARS-UCD1.2 and UMD3.1 reference genomes using a workflow based upon https://github.com/SichongP/SRA_variant_ search. The NCBI Genome Remapping Service (https://www.ncbi.nlm.nih.gov/genome/tools/remap) was used to translate variant positions between genome builds. The pileup function from sra-tools (https://github.com/ncbi/sra-tools, version 2.9.1) was used to generate pileup files at candidate variant loci. Subsequent genotyping was carried out for individuals with at least 10 reads per locus at the candidate loci.
For verification of the leading candidate SNP, 762 additional Hereford samples (semen, EDTA blood, or hair) were obtained. Of these samples, 289 were banked DNA samples from the American Hereford Association stored at Neogen GeneSeek (Lincoln, NE, USA), which were received in the form of isolated DNA. DNA was isolated from the blood and semen samples as described above. DNA from hair was isolated as previously described [12].
One of three genotyping methods was employed for validation of the candidate variant. Kompetitive Allele Specific PCR (KASP) genotyping was conducted using primers and probes designed with the KASP on Demand utility (LGC Genomics, Teddington, Middlesex, UK; Table S1). All KASP reactions were performed in duplicate on a CFX384 Touch Real-Time PCR machine (Bio-Rad Laboratories, Hercules, CA, USA) following the manufacturer's protocol. Non-template (negative) controls, three homozygous reference controls, three heterozygous controls, and four homozygous variant controls were run on each plate. The fluorophores HEX and FAM labeled the wildtype and variant probes, respectively. Results were visualized in CFX Maestro Software (Bio-Rad Laboratories, Hercules, CA, USA).
Any sample failing to genotype in duplicate via KASP was genotyped by either Sanger sequencing or droplet digital PCR (ddPCR). For ddPCR assays, primers to amplify a 136-bp fragment containing the candidate variant were designed in Primer3 [13] (Table S1). PrimeTime, double-quenched (ZEN/IowaBlack FQ) competitive probes were constructed to distinguish between the wildtype (T) and variant (C) alleles, labeled with HEX and FAM fluorophores, respectively. ddPCR reactions were performed using standard protocol on a QX200 ddPCR system (Bio-Rad Laboratories, Hercules, CA, USA). A non-template (negative) control was included as well as positive controls, which consisted of DNA from two affected calves, two obligate carriers, and two unaffected individuals. Samples of interest for allele quantification were run in duplicate or triplicate. The reaction included 1X ddPCR Supermix for Probes (no dUTP; Bio-Rad Laboratories, Hercules, CA, USA), 25 ng (10 ng/µL) of template DNA, primers at 0.18 µM each, probes at 0.02 µM each, and molecular grade water to a final reaction volume of 22 µL. Reactions were converted into approximately 20,000 one-nanoliter droplets using the QX200 Droplet Generator. Thermocycling included 10 min of enzyme activation at 95 • C and 39 cycles of denaturation (94 • C, 30 s) followed by annealing/extension (64 • C, 1 min). Enzyme deactivation (98 • C, 10 min) concluded the cycle. The plate was read in the QX100 Droplet Reader (Bio-Rad Laboratories, Hercules, CA, USA) and results were analyzed using QuantaSoft Software (Bio-Rad Laboratories, Hercules, CA, USA). Power to distinguish alleles was determined from the false-negative rate of the controls (e.g., power for detection of a variant allele = 1 − (wildtype droplets/total positive droplets), when template representing only the variant allele was provided).
Sanger sequencing was performed for 6 different variants at ACGT, Inc. (Wheeling, IL, USA) after PCR with primers designed in Primer3 [13] (Table S1) (Table S1) for 30 s, 72 • C for 45 s, a final extension at 72 • C for 10 min, then a 10 • C hold. PCR product cleanup was performed using 0.75 µL ExoSAP-IT (Applied Biosystems, Foster City, CA, USA) per 4 µL PCR product. Thermal cycling conditions consisted of a cycle at 37 • C for 30 min, 80 • C for 15 min, and 15 • C hold.
The cytochrome CYP26C1 (Bos taurus) structural model was generated based on the X-ray crystal structure of cyanobacterial CYP120A1 [14] (PDB: 2VE3) using sequence homology modeling program Modeller [15]; their protein sequences have a 31.4% identify and 46.8% similarity. The secondary structural predictions were conducted using Jnet online with the window size setup for 22 amino acids based on the structure of the α-helix segment from 179 to 200 amino acids. Conservation of the amino acid altered by the candidate SNP was also evaluated across multiple species using Multialign (http://multalin.toulouse.inra.fr/multalin).
Whole-genome sequence data generated for this project are available on the NCBI SRA (BioProject ID: PRJNA663547). Sanger sequencing data representing each CYP26C1 genotype have been deposited in GenBank (Accessions: MW123048 and MW123049). Novel candidate variants have been deposited in the European Variation Archive (EVA) (Project ID: PRJEB40605).

Pathologic Characteristics of Affected Calves
At autopsy, the 5 MD calves weighed between 32 and 41 kg. Unique and consistent hallmarks of the condition included bilateral skin tags 2-10 cm caudal to the commissure of the lips ( Figure 1A,B). The tags were 0.5-2.0 cm long with a central dermal core that attached through a short, 1 cm dermal band to cartilage ( Figure 1C). This cartilage extended to its origin at the zygomatic process of the temporal bone ( Figure 1D). The cartilage, from its origin, extended cranially for 3-5 cm and was encased in bone. The bony processes were 1.5-2.0 cm wide, 1.0 cm thick, and separated from underlying bones of the face by a 0.5 cm gap. A short, 1.0 cm dorsal lateral protrusion of bone at the origin of the bony process was also present. The bone-wrapped Meckel's cartilage was bilateral and consistent in each affected calf. When the bone was fractured during autopsy, a perfectly round 0.35 mm diameter tube of cartilage freely separated from the center of the bone ( Figure 1E). Several affected calves had additional skin tags near or several centimeters below the external acoustic meatus not associated with cartilage or bone ( Figure 1A,B). Additional coexisting facial deformities included megastomia in three, camplyognathia of the mandible in two, campylognathia involving mandible and maxilla in one, and brachgnathia inferior in three calves (Table 1). A single calf had a cleft palate. The calves with maxillary campylognathia also had asymmetry of the orbits with one located approximately 3 cm caudal to the other. All calves had hypoplasia of the masseter and temporalis muscles and pinnae that were low set and drooped.
Genes 2020, 11, x FOR PEER REVIEW 5 of 14 located approximately 3 cm caudal to the other. All calves had hypoplasia of the masseter and temporalis muscles and pinnae that were low set and drooped.
.   The dash (-) indicates an attribute that was not examined.
Histologic evaluation of the bony process originating from the temporal bone revealed a cartilage (Meckel's) core sandwiched between plates of bone ( Figure 1F). The bone and cartilage were surrounded by thick fibrous periosteum. No evidence of endochondral ossification of this cartilage was noted. Sections of eye, kidney, liver, brain, adrenal, spleen, skeletal muscle, thymus, intestine, and lymph node appeared normal. The skin tag had a thick dermal core with redundant collagenous stroma and most sections included normal adnexa and sinus hairs.
After disclosure of the MD defect to the breed association membership, three additional affected calves were reported, with the phenotype of two confirmed via digital image evaluation, and one by evaluation at the UNL-VDC. All affected calves had the presumed founder in their pedigree (Figure 2). The dash (-) indicates an attribute that was not examined.
Histologic evaluation of the bony process originating from the temporal bone revealed a cartilage (Meckel's) core sandwiched between plates of bone ( Figure 1F). The bone and cartilage were surrounded by thick fibrous periosteum. No evidence of endochondral ossification of this cartilage was noted. Sections of eye, kidney, liver, brain, adrenal, spleen, skeletal muscle, thymus, intestine, and lymph node appeared normal. The skin tag had a thick dermal core with redundant collagenous stroma and most sections included normal adnexa and sinus hairs.
After disclosure of the MD defect to the breed association membership, three additional affected calves were reported, with the phenotype of two confirmed via digital image evaluation, and one by evaluation at the UNL-VDC. All affected calves had the presumed founder in their pedigree ( Figure 2).

Whole-Genome Sequencing
Whole-genome sequence data averaged 13.0X coverage (range 10.8 to 15.5; standard deviation 1.32) and included 20 individuals (3 affected calves, parents of the sentinel cases (N = 6), parents of other affected calves (N = 4; including that of the 2018 calf) and 7 related individuals) ( Figure S1). 143 variants matched the hypothesized recessive mode of inheritance; all but 5 had RefSNP (rs) identifiers. 134 of the candidate variants were located between 10.3 and 15.9 Mb on chromosome 26 (Table S2). Evaluating the variants homozygous in the three affected calves without regard to the genotype of other individuals identified an extended region of homozygosity that also included this region of chromosome 26 (Table S3).
Additional evaluation was performed using WGS data from 101 animals sequenced for other projects in the lab (including 6 purebred Hereford and 8 Hereford-crosses), from 1577 cattle sequenced previously as part of the 1000 Bull Genomes Project [16], and 128 other Herefords provided by the American Hereford Association. These data resulted in the elimination of 141 of the 149 candidate variants as one or more of the additional animals were homozygous for the alternative allele. Of the remaining eight variants, one was a missense variant in Exon 3 of CYP26C1 (Chr26 g. 14404993T>C; ss7148511443), one SNP was located in the 5 untranslated region of TBC1D12 (Chr26 g. 15898152C>T; ss7148511444), and the other 5 were intronic or intergenic ( Table 2). The missense variant was the only variant with a predicted impact (SIFT score = 0, deleterious).

Candidate Variant Filtering
A query of the SRA yielded results for 3191 cattle at one or more of the 8 variants queried. These samples represented Hereford, Angus, Red Angus, and Simmental, among other breeds. At 6 of the remaining 8 variants queried, all cattle in the SRA search were homozygous reference, providing no information from which to further filter them from the data set. Two variants were eliminated as possible causative variants due to their presence in other breeds as MD appears to be unique to the Hereford and associated with one specific bloodline. For the first, Chr26 g. 10713132G>A, 6 of 1018 samples with genotype data were heterozygous for the alternative allele (frequency: 0.006). These 6 samples were of breeds including Romagnola, Simmental, and Original Braunvieh. At Chr26 g. 10794674G>A, 7 of 1083 samples were heterozygous and 1 homozygous for the alternative allele (frequency: 0.0083). These 8 samples were also Romagnola, Simmental, and Original Braunvieh.
Sanger sequencing was performed on the 6 remaining candidate variants in one of the additional affected calves, three additional parents of affected calves, the presumed founder bull, the sire of the presumed founder, and 25 animals with the presumed founder in both their maternal and paternal pedigree. These data ruled out 4 of the 6 remaining candidate variants due to identification of animals that genotyped homozygous for the alternative allele but did not display the MD phenotype or due to the MD affected calf not genotyping as homozygous for the alternative allele. Of the two variants that could not be ruled out with additional data (Chr26 g. 14404993T>C and Chr26 g. 15898152C>T), due to its predicted function, the missense variant in CYP26C1 was the primary candidate of interest and focused upon for further analyses.

CYP26C1 Variant Genotyping
782 Hereford cattle were genotyped for the CYP26C1 SNP including the WGS animals (N = 20), other cattle from two of the reporting herds (N = 623), the suspect founder, two additional affected calves, and parents of those affected calves (N = 4) ( Table 3). Of these 782, 624 were homozygous for the reference allele, 153 heterozygous, and the 5 affected calves homozygous for the variant allele. Within these Herefords genotyped, 106 had the presumed founder bull in both the sire and dam side of their pedigree, 327 had the presumed founder bull in one side of their pedigree, and 348 animals were not descendants of the presumed founder (Table 3). Notably, the suspect founder genotyped heterozygous, all additional affected calves from which DNA was available were homozygous for the alternative allele, and parents of the additional affected calves were heterozygous (DNA on one dam not available). Twelve animals were calves of both a sire and dam that were confirmed carriers of the CYP26C1 variant allele; of these, 2 were homozygous for the reference allele, 5 heterozygous, and 5 homozygous for the variant allele (all affected with the described MD phenotype). In total, 3371 genotypes for the CYP26C1 g. 14404993T>C locus were examined (Table 4). With the exception of 17 of the 20 cattle sequenced for this study, the variant was not identified in any WGS-derived genotype.
Power to distinguish between wildtype and variant alleles by ddPCR was 0.999 and 0.998, respectively. An average of 9999.9 droplets (st dev = 1599.8) were read per sample, with each run in 2 or more replicates. The germ line variant allele frequency of the maternal grandsire and great grandsire of the suspect founder did not exceed that of wildtype controls (Table S4), refuting a hypothesized mosaic origin of the variant in either bull.  Table 3.

Predicted Impact on Protein Function
Alignment of the amino acid sequence of CYP26C1 indicated that the residue altered (p. L188P) was conserved across B. taurus, Ovis aries, Sus scrofa, Equus caballus, Homo sapiens, Mus musculus, Danio rerio, Xenopus tropicalis, and Canis lupus familiaris (Figure 3). Structural modeling demonstrated the L188P mutation site is located in an α-helix segment from 179 to 200 amino acids. This is not within the heme active site at which the substrate, all-trans-retinoic acid (atRA), binds. Secondary structural predictions, however, indicate the α-helix will be broken if leucine (L) is substituted with proline (P) at position 188 ( Figure 4).
Genes 2020, 11, x FOR PEER REVIEW 9 of 14 Power to distinguish between wildtype and variant alleles by ddPCR was 0.999 and 0.998, respectively. An average of 9999.9 droplets (st dev = 1599.8) were read per sample, with each run in 2 or more replicates. The germ line variant allele frequency of the maternal grandsire and great grandsire of the suspect founder did not exceed that of wildtype controls (Table S4), refuting a hypothesized mosaic origin of the variant in either bull.

Predicted Impact on Protein Function
Alignment of the amino acid sequence of CYP26C1 indicated that the residue altered (p. L188P) was conserved across B. taurus, Ovis aries, Sus scrofa, Equus caballus, Homo sapiens, Mus musculus, Danio rerio, Xenopus tropicalis, and Canis lupus familiaris (Figure 3). Structural modeling demonstrated the L188P mutation site is located in an α-helix segment from 179 to 200 amino acids. This is not within the heme active site at which the substrate, all-trans-retinoic acid (atRA), binds. Secondary structural predictions, however, indicate the α-helix will be broken if leucine (L) is substituted with proline (P) at position 188 ( Figure 4).

Discussion
In this investigation, MD of Hereford calves was associated with a missense mutation in CYP26C1. MD is detrimental to the viability of affected animals due to an impaired ability to nurse. The identification of this associated and putatively causative marker provides breeders the ability to test for carriers to avoid mating decisions that may produce an affected calf. Animals with the presumed founder on both sides of their pedigree originated from over 40 different breeders showing the widespread implications of MD within the Hereford breed. Further, the implication of CYP26C1 provides insight into the role of this gene in the regulation of RA during the development of craniofacial structures.
Facial structures including the mandible, maxillary musculature, and some structures of the inner ear and hindbrain are derived embryonically from CNCC that migrate from the neural tube to PA1. This region of the embryo is the most rostral of the anteroposteriorally segmented structures present in early development [17]. Expression of Dlx homeobox genes in the CNCC of PA1 provides patterning information that results in specification of the mandibular and maxillary arches [4,6]. Differentiation between the mandible and maxilla is driven by ET1 signaling in PA1 that activates Dlx5 and Dlx6 resulting in formation of Meckel's cartilage in the mandibular, but not the maxillary prominence [5,6]; Meckel's cartilage, a transient cartilaginous plate, serves as a scaffold during formation and elongation of the intramembranous mandible bone [18].
Migration and proliferation of CNCC is guided by gradients of signals including fibroblast growth factor, retinoic acid, and Wnt [6,19]. Retinoic acid (RA) plays a critical role in PA1 development because it acts as a repressor of ET1 [5]. Excess RA disrupts ET1 signaling acting on the Dlx genes. This excess exposure of PA1 to RA results in teratogenesis of the mandible, specifically altering development of the Meckel's cartilage [5]. Similar effects were observed in both Dlx5/Dlx6

Discussion
In this investigation, MD of Hereford calves was associated with a missense mutation in CYP26C1. MD is detrimental to the viability of affected animals due to an impaired ability to nurse. The identification of this associated and putatively causative marker provides breeders the ability to test for carriers to avoid mating decisions that may produce an affected calf. Animals with the presumed founder on both sides of their pedigree originated from over 40 different breeders showing the widespread implications of MD within the Hereford breed. Further, the implication of CYP26C1 provides insight into the role of this gene in the regulation of RA during the development of craniofacial structures.
Facial structures including the mandible, maxillary musculature, and some structures of the inner ear and hindbrain are derived embryonically from CNCC that migrate from the neural tube to PA1. This region of the embryo is the most rostral of the anteroposteriorally segmented structures present in early development [17]. Expression of Dlx homeobox genes in the CNCC of PA1 provides patterning information that results in specification of the mandibular and maxillary arches [4,6]. Differentiation between the mandible and maxilla is driven by ET1 signaling in PA1 that activates Dlx5 and Dlx6 resulting in formation of Meckel's cartilage in the mandibular, but not the maxillary prominence [5,6]; Meckel's cartilage, a transient cartilaginous plate, serves as a scaffold during formation and elongation of the intramembranous mandible bone [18].
Migration and proliferation of CNCC is guided by gradients of signals including fibroblast growth factor, retinoic acid, and Wnt [6,19]. Retinoic acid (RA) plays a critical role in PA1 development because it acts as a repressor of ET1 [5]. Excess RA disrupts ET1 signaling acting on the Dlx genes. This excess exposure of PA1 to RA results in teratogenesis of the mandible, specifically altering development of the Meckel's cartilage [5]. Similar effects were observed in both Dlx5/Dlx6 and ET1 receptor mutant mice [20], supporting the importance of RA in the development of mandibular structures.
Dysregulation of the RA gradient has the most severe effect on craniofacial morphology at the 9 to 14 somite stage of development [5]; during this period, migrating CNCC reach PA1 [21]. Nine-somite mouse embryos treated with excess RA had malformed Meckel's cartilage [4,5] while the cartilage developed normally when treated outside of this developmental time point [5]. For cattle, this critical developmental window would fall around day 21 or 22 of gestation [22,23].
The concentration of RA in the developing embryo is regulated in a tissue-specific manner by opposing actions of synthesizing (Raldh) and metabolizing (Cyp26) enzymes [24]. The CYP26C1 gene and its counterparts, CYP26A1 and CYP26B1, oxidize RA into inactive polar metabolites enabling the maintenance of the RA gradient across the pharyngeal arches [25]. Studies elucidating the expression patterns of CYP26C1 show it is expressed in PA1 during the critical time point for RA teratogenesis [24].
There are multiple lines of evidence that alteration of CYP26C1 function results in developmental abnormalities. Knockdown of CYP26C1 in zebrafish embryos resulted in increased RA levels associated with atypical development of the pharyngeal arches and otic vesicles [26]. CYP26A1 and CYP26C1 null mice also had deficient development of PA1 and PA2 and altered migration of CNCC [27]. These studies indicate that normal function of CYP26C1 is critical to avoid morphological changes in craniofacial structures.
We postulate that the association of the CYP26C1 missense mutation with abnormal development of the jaw in MD calves can be attributed to the loss of catalytic activity resulting in reduced ability of the enzyme to metabolize RA. The affected amino acid is not located in the canonical heme thiolate-binding motif (FxxGxxxCxG) that is characteristic of the P450 family of proteins to which CYP26C1 belongs [28]. However, Kim and Kang [29] showed proline is a helix breaker and our model supports that the L188P mutation results in a structural change of the α-helix. This change may affect the active site structure for substrate binding and cause a subsequent loss of function impairing metabolism of RA. Further, this amino acid residue is conserved across species. Both the protein model and the evolutionary maintenance of the amino acid serve as evidence the CYP26C1 L188P mutation is deleterious to protein function.
Reports of craniofacial deformities in cattle commonly include cleft palate or more severe dysplasias [30][31][32]. A recent report from India identified a calf with a pathology similar to MD [33]. This, however, is the first comprehensive report of the novel MD phenotype in Herefords, which we believe is attributed to a de novo mutation in the sire common to the paternal and maternal pedigrees of all affected calves. Semen from the suspect founder genotyped heterozygous for the CYP26C1 SNP with no evidence of germ-line mosaicism as studied by ddPCR. Semen of this bulls' sire and both maternal grandsires were homozygous for the reference allele at the CYP26C1 SNP, also confirmed by ddPCR. DNA was not available from the (deceased) mother of the presumed founder bull, therefore origination of the variant in the maternal lineage of the presumed founder cannot be excluded. Even with the possibility of maternal origin, the bull common to all pedigrees of affected calves has been the primary perpetuator of the variant resulting in the manifestation of this phenotype. To this point, no carrier animals have been identified that do not have the presumed founder in their pedigree. Further, both the presumed founder and his sire genotyped as heterozygous for the only other remaining candidate variant (Chr26 g. 15898152C>T). If that variant was instead causative, we would expect to have identified MD cases in descendants of the sire of the suspect founder without a pedigree connection to the suspect founder himself; that has not been the case, contributing to the body of evidence that the variant in CYP26C1 is causative of this condition.
Several human conditions exhibit similar pathologies to MD calves. One hallmark of MD in all affected calves is the presence of skin tags located along the fusion site of the maxillary and mandibular prominences. Strikingly similar skin tags are observed in Hemifacial Microsomia patients [3]. Human focal facial dermal dysplasia, Type IV, is also characterized by skin lesions and polyps on the buccal mucosa located at the same fusion site [28]. Notably, both human syndromes are attributed to similar developmental pathways as implicated in MD. Focal facial dermal dysplasia is caused by loss of function mutations in CYP26C1 [28] and one model of causation for hemifacial microsomia includes altered migration, differentiation, and proliferation of CNCC in PA1 [3].
A condition also termed mandibulofacial dysostosis is described in human literature as a heterogeneous anatomic group of disorders; the described condition in Hereford cattle is not analogous to human MD. The human MD condition more broadly affects chondrocyte and osteoblast differentiation with dysostosis most apparent in, but not limited to, the face [34]. Implicated mutations in EFTUD2 in human MD are unrelated to CNCC migration and differentiation [35,36]. Further, the normal facial features of humans and bovines differ dramatically with cattle being more dolichocephalic; thus, we do not anticipate homologous anatomic outcomes due to disruptions in branchial arch related processes.
A defining and novel characteristic of Hereford MD is the retention of Meckel's cartilage past early development and shortening of the mandible. No human conditions have been described exhibiting this precise pathology, but similarities in a group of human syndromes termed "retinoic acid embryopathies" exist. These embryopathies have been described to result in micrognathia, cleft palate, and microtia or anotia [37]. One of the affected MD calves exhibited cleft palate as observed in human retinoic acid embryopathies. Additionally, hemifacial microsomia patients manifest hypoplasia of the mandible due to decreased blood supply to the Meckel's cartilage [38].
Finally, calves exhibited a range of severity in phenotype in the truncation of the lower jaw, presence of additional skin tags, and in one case the additional feature of a cleft palate. Variation in human syndromes affecting craniofacial development is also common. For example, the severity of hemifacial microsomia can range from moderate hypoplasia to complete absence of a portion of the jaw with other neural and muscular symptoms [4]. Similarly, in a set of infants with malformations due to RA exposure, 14% exhibited cleft palate [37]. We propose that variation in phenotype of the MD calves may be attributed to intrinsic variability in RA availability during development.