Identification of Inherited Retinal Disease-Associated Genetic Variants in 11 Candidate Genes

<jats:p>Inherited retinal diseases (IRDs) display an enormous genetic heterogeneity. Whole exome sequencing (WES) recently identified genes that were mutated in a small proportion of IRD cases. Consequently, finding a second case or family carrying pathogenic variants in the same candidate gene often is challenging. In this study, we searched for novel candidate IRD gene-associated variants in isolated IRD families, assessed their causality, and searched for novel genotype-phenotype correlations. Whole exome sequencing was performed in 11 probands affected with IRDs. Homozygosity mapping data was available for five cases. Variants with minor allele frequencies ≤ 0.5% in public databases were selected as candidate disease-causing variants. These variants were ranked based on their: (a) presence in a gene that was previously implicated in IRD; (b) minor allele frequency in the Exome Aggregation Consortium database (ExAC); (c) in silico pathogenicity assessment using the combined annotation dependent depletion (CADD) score; and (d) interaction of the corresponding protein with known IRD-associated proteins. Twelve unique variants were found in 11 different genes in 11 IRD probands. Novel autosomal recessive and dominant inheritance patterns were found for variants in Small Nuclear Ribonucleoprotein U5 Subunit 200 (SNRNP200) and Zinc Finger Protein 513 (ZNF513), respectively. Using our pathogenicity assessment, a variant in DEAH-Box Helicase 32 (DHX32) was the top ranked novel candidate gene to be associated with IRDs, followed by eight medium and lower ranked candidate genes. The identification of candidate disease-associated sequence variants in 11 single families underscores the notion that the previously identified IRD-associated genes collectively carry &gt; 90% of the defects implicated in IRDs. To identify multiple patients or families with variants in the same gene and thereby provide extra proof for pathogenicity, worldwide data sharing is needed.</jats:p>

section_content: Inherited retinal diseases (IRDs) comprise of a clinically and genetically heterogeneous group of hereditary retinal degenerations, which affect around 1 in 3000 people. Inherited Retinal Diseases may occur both in syndromic and non-syndromic forms, and can be inherited in all the Mendelian inheritance patterns. Inherited Retinal Diseases are one of the leading causes of vision loss in young individuals in the developed countries, with a significant impact on patients and society [1] [2] [3]. In the past two decades, new therapies have been developed that rely on knowledge of the underlying genetic defects [4] [5] [6]. 

The high genetic diversity and overlapping phenotypes of IRDs pose significant challenges in the genetic elucidation of these disorders. The spectrum of genotyping technologies employed in the last three decades, ranging from positional cloning, (mouse) candidate disease gene analyses, linkage analysis, homozygosity mapping, and whole exome sequencing (WES) have led to a wide body of knowledge regarding the molecular genetics of IRDs. Currently, mutations in 261 genes have been discovered to be involved in IRDs [7]. Nevertheless, there is a trend in which the more recently discovered IRD-associated genes are responsible for a very small percentage of cases, in comparison to the genes that were discovered in the earlier studies. Two WES studies in the Saudi-Arabian and the Dutch population showed that only ~3-4% of IRD cases carry variants in novel IRD-associated genes [8] [9] [10] [11] [12]. Consequently, finding additional cases carrying variant(s) in the same candidate gene is often challenging. 

A large cohort of IRD cases is crucial to assess the contribution of unique variants in causing diseases. The European Retinal Diseases Consortium (ERDC) [13], consisting of 17 research groups, was formed as a collaborative effort to elucidate molecular genetics causes of IRDs. Within the ERDC, information on candidate IRD genes was shared which enabled the identification of multiple cases or families with variants in the same gene [14] [15] [16] [17]. 

In this study, we report on the identification of potential candidate disease variants in 11 genes that were found in single families. We assessed the causality of these candidate disease variants, searched for novel genotype-phenotype correlations and investigated the occurrence of different inheritance patterns for variants in known IRD-associated genes.
section_name: Introduction

section_name: Materials and Methods

section_content: The study protocol adhered to the tenets of the declaration of Helsinki and received approval from the respective local ethics committees (2011, MEC-2010-359, ABR-NL34152. 078. 10; 2004, 03/362 17/YH/0365; 2015, 176953, 15/YH/0365; 2008-A01238-47). Written informed consent was obtained from all participants or parents of children prior to their inclusion in this study. The 11 probands belong to eight autosomal recessive retinitis pigmentosa (RP) families, two autosomal dominant RP families, one autosomal recessive cone-rod dystrophy (CRD) family, one achromatopsia (Figures 1 and 2, Tables 1 and 2 ). Each patient was reviewed for clinical examination including visual acuity, fundus photography, fundus autofluorescent, and optical coherence tomography.
section_name: Subject and Clinical Evaluation

section_content: Genomic DNA was isolated from peripheral blood lymphocytes according to standard salting out procedures. Whole exome sequencing analysis consisting of targeted capture, sequencing platform, and mapping was independently performed by each center using several types of software (Supplementary Table S1 ) [18] [19] [20]. Sequences were aligned to National Center for Biotechnology Information (NCBI) Build 37 of the human genome [21] using corresponding mapping tools and copy number variant analyses were performed. Single nucleotide polymorphism (SNP) microarray analysis was performed in five families and the runs of homozygosity were assessed using PLINK or Homozygosity mapper [20, 22].
section_name: Genetic Analysis

section_content: Whole exome sequencing identified approximately 120,000 variants in each proband. Inclusion of rare variants with minor allele frequency (MAF) ≤ 0. 5% was performed as an initial filtering step. Thereafter, protein coding and splice site variants were filtered and assessed using in silico prediction tools to prioritize the remaining variants. Four parameters were used to assess the pathogenicity of candidate variants: (a) minor allele frequencies obtained from the Exome Aggregation Consortium (ExAC) [23] database containing 60,706 exomes of unrelated healthy individuals, accessed September 2017, (b) in silico pathogenicity assessment using the combined annotation dependent depletion (CADD) score, (c) interaction with known IRD-associated proteins, and (d) known association of the relevant gene with IRDs. Using these criteria, all of the candidate genes were scored and ranked according to the number of matching parameters. The total score corresponded to the likelihood of a gene to be associated with IRDs. Yet Another Scientific Artificial Reality Application (Yasara) [24] was used to model the possible three-dimensional (3D) structural alterations caused by Small Nuclear Ribonucleoprotein U5 Subunit 200 (SNRNP200) amino acid changes [25]. Primers for candidate variants confirmation were designed with Primer3 plus [26].
section_name: Variant Prioritization

section_content: Each ERDC partner shared candidate genes in biannual meetings. Thereafter, the corresponding genes were screened for the presence of variants in WES data and homozygosity mapping data were assessed for homozygous regions larger than 4 Mb containing these candidate genes. In addition,
section_name: Data Sharing

section_content: Whole exome sequencing identified approximately 120,000 variants in each proband. Inclusion of rare variants with minor allele frequency (MAF) ≤ 0. 5% was performed as an initial filtering step. Thereafter, protein coding and splice site variants were filtered and assessed using in silico prediction tools to prioritize the remaining variants. Four parameters were used to assess the pathogenicity of candidate variants: (a) minor allele frequencies obtained from the Exome Aggregation Consortium (ExAC) [23] database containing 60,706 exomes of unrelated healthy individuals, accessed September 2017; (b) in silico pathogenicity assessment using the combined annotation dependent depletion (CADD) score; (c) interaction with known IRD-associated proteins; and (d) known association of the relevant gene with IRDs. Using these criteria, all of the candidate genes were scored and ranked according to the number of matching parameters. The total score corresponded to the likelihood of a gene to be associated with IRDs. Yet Another Scientific Artificial Reality Application (Yasara) [24] was used to model the possible three-dimensional (3D) structural alterations caused by Small Nuclear Ribonucleoprotein U5 Subunit 200 (SNRNP200) amino acid changes [25]. Primers for candidate variants confirmation were designed with Primer3 plus [26].
section_name: Variant Prioritization

section_content: Each ERDC partner shared candidate genes in biannual meetings. Thereafter, the corresponding genes were screened for the presence of variants in WES data and homozygosity mapping data were assessed for homozygous regions larger than 4 Mb containing these candidate genes. In addition, candidate genes were posted at the open access ERDC website to increase the probability of finding additional cases carrying variants in the same candidate genes [13].
section_name: Data Sharing

section_content: We performed an assessment to predict the pathogenicity of selected candidate variants by combining: (a) minor allele frequency (ExAC); (b) in silico prediction program (CADD score); (c) whether there is any interaction with retinal-associated protein (from previous studies or STRING protein network database) [27]; and, (d) whether there is any described association with IRDs. Subsequently, the total score was used to consider the likelihood of variants to be causative. In addition, we also implemented the American College of Medical Genetics (ACMG) guidelines to predict variants based on their pathogenicity [28]. Upon this assessment, variants in the known IRD-associated gene SNRNP200 is the highest ranked classification 'score 4'. For the other genes, two genes have a 'score 3' using our assessment, i. e., DEAH-Box Helicase 32 (DHX32) and Zinc Finger Protein 513 (ZNF513). The eight remaining genes have a lower pathogenicity scoring with 'score 2' or low pathogenicity scoring of 'score 1'. Hence, these genes are better candidate in comparison to the other genes. In addition, we provide the amino acid conservation for all the missense variants identified in this study (Figure 3 ). Clinical details of affected individuals of all families are described in Table 2.
section_name: Results

section_content: In the consanguineous Pakistani family G, we performed microarray SNP analysis prior to WES. In the largest homozygous region (14. 6 Mb), we identified a homozygous c. 3269G>A (p. (Arg1090Gln)) variant in SNRNP200. This variant occurs in 3 of 121,400 (0. 000025) alleles in ExAC and is absent in 180 ethnically matched in-house controls. The substitution affects a highly conserved amino acid and has a PhyloP score [29] of 5. 77 (range-14. 1; 6. 4). The new residue has a small physicochemical difference in comparison to the wild type amino acid, with a Grantham [30] score of 43 (range 0-215). The wild type arginine residue in position 1090 is located in the ratchet helix of the helical bundle (HB) of human Bad Response to Refrigeration 2 (hBRR2) that forms a unidirectional tunnel for RNA molecules to promote RNA unwinding events (Figure 4A, B ). Arg1090 seems to be in direct contact with the RNA molecule and an alteration in this position may disrupt the unwinding of hBRR2.
section_name: Score 4 Candidate IRD Gene

section_content: The DHX32 extension variant c. 2176del (p. (Glu726Asnfs*57)) in family A will result in a frameshift with read-through of the stop codon resulting in a transcript encoding 57 additional amino acids when compared to the wild type protein. This variant has a minor allele frequency of 0. 000016 in the ExAC database. DHX32 is one of DEAD (Asp-Glu-Ala-Asp) box proteins with a putative role as a RNA helicase that is involved in nuclear regulation and mitochondrial gene expression. As the variant occurs in the last exon of DHX32 there will be no nonsense mediated decay. The addition of 57 amino acids could however destabilize the protein's folding or it might interfere in protein-protein interactions. 

The third candidate gene is ZNF513. The variant c. 724C>T (p. (Arg242Cys)) was identified in two individuals from family K affected with late onset autosomal dominant RP. This variant is absent in the ExAC database. The variant affects a highly conserved amino acid and substitutes a nucleotide with a PhyloP score of 2. 71. In addition, the wild type arginine residue has a large physicochemical difference with cysteine, with a Grantham score of 180. Based on the family history (grandfather deceased at age 44 years when disease may not yet have appeared; great-grandfather mentioned to be blind) there is strong supporting evidence for a dominant inheritance pattern. To further rule out autosomal recessive inheritance, a detailed analysis of copy number variants was performed. In the fully covered gene, there was no evidence for the presence of deletions or duplications in ZNF513 in the affected proband.
section_name: Score 3 Candidate IRD Genes

section_content: In family B, we identified a homozygous missense c. 925G>A (p. (Ala309Thr)) variant in DSCAML1, a gene that encodes the Down syndrome cell adhesion molecule-like 1 transmembrane protein. 

The substitution occurs in the amino acid that is highly conserved. The variant affects a highly conserved nucleotide with a PhyloP score of 5. 86, but the altered amino acid is expected to have a relatively small physicochemical difference with a Grantham score of 58. This variant is found in 0. 000022 of alleles in the ExAC database. All of the in silico prediction programs predicted this substitution to be deleterious, supported by the high CADD score of 34. The variant is located in the region that encodes the Immunoglobulin C2-type domain of the gene. Down syndrome cell adhesion molecule (DSCAM) and DSCAML1 are known to be cell adhesion molecules that play roles in several neurodevelopmental functions, including neuronal self-avoidance [33]. This adhesion molecule is expressed in a subset of amacrine cells and most retinal ganglion cells (RGCs) soon after they develop and promotes laminar specificity in the retinal synaptic connections via homophilic intracellular interactions [33]. A recent study in mice shows the expression of Dscaml1 in the rod circuit, and the lack of Dscaml1 leads to fasciculate rod bipolar cell dendrites and clumped amacrine cell bodies, suggesting its role in the neuronal self-avoidance. Moreover, this study also suggested the importance of Dscaml1 for rod bipolar cell dendrite organization and amacrine cell mosaic patterning [34]. 

Segregating compound heterozygous variants c. 2318C>T and c. 2833C>T (p. (Pro773Leu) and p. (Arg945*)) were identified in Otogelin Like (OTOGL) in family F with non-syndromic autosomal recessive RP. The missense variant alters a highly conserved nucleotide and amino acid (proline), with a conserved nucleotide with PhyloP score of 4. 00. This variant occurs in 89/118,882 alleles (0. 00075) in the ExAC database, i. e., heterozygously in 85 persons and homozygously in two persons. The allele frequency of the nonsense variant was 0. 000076 in the ExAC database and predicted to result in nonsense mediated decay. If there would be a shortened protein, it will no longer include two von Willebrandt factor, type D domains, two C8-domains, and the C-terminal domain. Pathogenic mutations in OTOGL are associated with sensorineural hearing loss [35, 36]. Based on their allele frequencies, we considered both variant to be less likely associated with the disease phenotype; however, no other outstanding variants were observed. Homozygosity mapping with the 250. 000 SNP array had indicated two homozygous regions in which no currently known IRD-associated genes appear. 

In the four-generation family J, we detected a heterozygous missense c. 1069G>A (p. (Ala357Thr)) variant in the E1A Binding Protein P300 (EP300) gene that inherited in an autosomal dominant manner. The variant alters a moderately conserved amino acid and affects a highly conserved nucleotide (PhyloP: 6. 34). This substitution is absent from the ExAC database and all in silico predictions show this variant to be pathogenic. It occurs in exon 4 that encodes the transcription adaptor putative Zinc finger domain of the EP300 protein. The EP300 gene encodes p300 histone acetyltransferase that modulate transcription through chromatin remodeling and also plays role in cell proliferation and differentiation [37]. EP300 is ubiquitously expressed in adult tissue. Heterozygous truncating variants in this gene have previously been associated with the Rubinstein-Taybi syndrome, a condition that is characterized by short stature, intellectual disability, facial dysmorphology, and often retinal disorders [38]. The non-ocular features of Rubinstein-Taybi syndrome were not observed in the proband. A study using conditional KO Ep300 mice in differentiating rods or cones using opsin driven Cre-recombinase showed that p300 is required for photoreceptor-specific chromatin organization, function, gene expression, and to maintain cell identity [39]. In addition, a recent study in zebrafish light-induced larvae suggested the role of EP300 in protecting photoreceptor cells from light-induced damage and shows that EP300 activation as a potential approach for the treating retinal degenerative diseases [40]. Taken together, this evidence suggests that EP300 is as plausible candidate gene that is associated with retinal phenotype in this family. Loss of function may lead to a severe phenotype, while missense variant in this family may still retain all the essential domain of p300, thus only causes retinal phenotype. 

Homozygosity mapping in combination with WES analysis detected a homozygous missense c. 930C>G (p. (Phe310Leu) variant in Farnesyl-Diphosphate Farnesyltransferase 1 (FDFT1) in two siblings from family C with rod-cone dystrophy. FDFT1 resides in the single homozygous region shared by the two affected siblings spanning 11. 6 Mb. The variant affects a highly conserved phenylalanine, and a poorly conserved nucleotide (PhyloP: 1. 09). The mutant leucine is also present in Baker's yeast (Supplemental Figure S1 ). This substitution occurs in 87 out of 121,400 alleles (0. 00072). The FDFT1 gene encodes farnesyl-diphosphate farnesyltransferase 1 enzyme, which catalyzes the dimerization of farnesyl-diphosphate to form squalene in the mevalonate pathway. Two IRD-associated genes have previously been linked to the process of prenylation [41]. The mevalonate kinase (MVK) gene has been implicated in non-syndromic RP [42], and mutations in Member RAS Oncogene Family 28 (RAB28), encoding a prenylated RAB28 protein, are a cause of CRD [17]. Despite this 'guilt by association', the information provided above suggest this variant to be benign. 

In family E with achromatopsia, we identified a homozygous c. 769C>T (p. (Tyr257His)) substitution in the NDRG Family Member 4 (NDRG4) gene. This gene resides in the largest overlapping homozygous region of 76. 5 Mb shared by the two affected siblings, but not by the healthy sibling. This highly conserved nucleotide variant (PhyloP: 4. 48) is absent in the ExAC database and alters a highly conserved tyrosine. The variant is predicted to be damaging by all prediction programs and occurs in the Ndr domain of the NDRG4. The gene is expressed in the retina based on the ocular tissue database, human proteome map, Gepis tissue and human retinal transcriptome. No association or interaction with retinal-associated protein is known for this protein. When compared to the missense variant identified in this study, the Ndrg4-deficient mice present with impaired spatial learning and memory, but normal motor function [43]. Previous studies have shown that Ndrg4 mRNA is expressed in the amacrine cells, horizontal cells, and ganglion cells in the developing mouse retina, and amacrine and possibly ganglion cells at 96 hours post fertilization in zebrafish [44]. 

A hemizygous c. 164G>T (p. (Ser55Ile)) variant in Synaptotagmin Like 4 (SYTL4) was found in family H, which was not present in the ExAC database. No other family members were available to investigate the X-linked mode of inheritance. The substitution alters an evolutionary highly conserved serine, and also affects a highly conserved nucleotide with a PhyloP score of 5. 21. The variant occurs in exon 4, which encodes the Rab binding domain. All the in silico tools suggest that the variant damages the protein. The larger size and the altered hydrophobicity of the new residue may disrupt the forming of the hydrogen bond between serine at position 55 with lysine at position 49 in the wild type situation. The altered residue is located on the protein surface, and thus may affect interaction with other molecules or protein domains. SYTL4 encodes a member of the synaptotagmin like protein family, which contains the N-terminal Rab27 binding domain and C-terminal tandem C2 domain. The encoded protein binds to Syntaxin Binding Protein 1 (STXBP1), Member RAS Oncogene Family 3A (RAB3A), Member RAS Oncogene Family 8A (RAB8A) and Member RAS Oncogene Family 27B (RAB27B), and interacts with Myosin VA (MYO5A). 

In family I, we identified a homozygous nonsense c. 418C>T (p. (Arg140*)) variant in Transmembrane P24 Trafficking Protein 7 (TMED7). This variant is absent from the ExAC database. TMED7 might play an important role in vesicular protein trafficking. This variant will result in either a truncated protein product lacking the EMP24_GP25L or GOLD domain, or may be subjected to nonsense mediated decay. 

Thus far, no disorder has been associated with mutations in this gene. Individuals carrying this variant demonstrated an early onset visual impairment in the first decade of life. In addition, nystagmus, developmental delay and tonic clonic seizures also occurred in the patient.
section_name: Score 2 Candidate IRD Genes

section_content: We identified a homozygous missense variant c. 4G>A (p. (Ala2Thr)) in G Protein-Coupled Receptor 45 (GPR45). This variant occurs in a homozygous state in two affected sisters with CRD in a consanguineous Israeli Muslim Arab family D. The variant is predicted to result in a small physicochemical change with a Grantham score of 58; the mutant threonine is larger than alanine, which may destabilize GPR45. The variant affects a moderately conserved nucleotide with a PhyloP score of 1. 90 and the alanine at position 2 is moderately conserved up to frog. This variant is not present in the ExAC database and all in silico tools predict the variant to be deleterious. In addition, the variant resides outside the region that encodes functional protein domain. GPR45 is a member of the G protein-coupled receptor (GPCR) superfamily that is known to modulate physiological processes in the neurological systems [45].
section_name: Score 1 Candidate Gene

section_content: Whole exome sequencing has been instrumental for IRD gene discoveries since 2011 [15, 46, 47]. The implementation of this technology has contributed to the identification of ~40 novel genes that are associated with syndromic and non-syndromic IRDs in the last six years. Furthermore, novel genotype-phenotype correlations were identified for several IRD-associated genes by using exome sequencing. A tendency can be observed, as the early-discovered genes are more frequently mutated than the more recently identified genes. Previous genetic studies employing homozygosity mapping and WES in the Saudi Arabian population shows a small proportion of IRD cases (3-4%) that carry a mutation in novel IRD-associated genes [8, 9]. Additionally, a gene panel study in 537 IRD cases demonstrated 79% of diagnostic yield in the genes that were discovered between 1995 and 2004 [48]. Together with challenges in interpreting novel variants, this low mutational load phenomenon led to the identification of many variants in novel candidate genes in singleton cases without any clear functional consequences. 

The ERDC aims to identify the remaining unidentified IRD-associated genes by data sharing. In the past few years, the ERDC was able to identify several genes associated with IRDs [14] [15] [16] 49]. To enable discovery of additional families harboring overlapping pathogenic variants the online Web-tool GeneMatcher may also be utilized [50]. This tool allows for researchers and clinicians to upload a gene with information about the identified genetic variants and the phenotype found in the patient that carries this variant or combination of variants. 

Upon screening of WES and identity-by-descent data, 11 genes remained as candidates that occur in single families or cases. Using our pathogenicity assessment SNRNP200, DHX32, and ZNF513 are the top ranked candidate genes based on the allele frequency, in silico prediction and interaction with the known IRD genes. Screening of IBD and WES data of unsolved IRD cases did not yield additional cases. The lower ranked candidate genes have generally weaker pathogenicity prediction scores, however all stand out when compared to other variants seen in the corresponding patients. 

Besides the identification of novel candidate genes, we also revealed a novel mode of inheritance in two previously described IRD-associated genes. The SNRNP200 variant is located in the Sec63 domain and predicted to be deleterious. Interestingly, an overlapping variant in this position, c. 3269G>T (p. (Arg1090Leu), Figure 4C ), has been previously identified in a large Chinese family segregating an autosomal dominant form of RP. SNRNP200 encodes BRR2, a member of the spliceosome complex. Functional assays in yeast using a Brr2p construct containing p. (Arg1090Leu) demonstrated a splice defect, suggesting the association of SNRNP200 with the disease phenotype. In the autosomal dominant RP associated variant, the wild type arginine is altered to a nonpolar and hydrophobic leucine. On the other hand, the autosomal recessive RP-associated variant alters to glutamine that is uncharged but has the same polarity as the wild type amino acid. The comparison of the two new residues in the 3D structure of SNRNP200 shows that the alteration to leucine shows a greater structural difference in comparison to glutamine. We hypothesize that p. (Arg1090Gln) only causes disease in a homozygous state, as its effect on the SNRNP200 activity is approximately half of the effect of p. (Arg1090Leu). The causality of autosomal recessive variants in SNRNP200 was suggested by Wang et al. In a Chinese Leber congenital amaurosis patient who carried a homozygous p. (Pro1045Thr) variant [51]. Moreover, Bujakowska and colleagues proposed autosomal recessive inheritance in a RP patient with compound heterozygous variants in SNRNP200 consisting of the same hypomorphic missense variant p. (Pro1045Thr) and a 1. 1-Mb deletion spanning SNRNP200 [52]. 

The ZNF513 mutant residue Cys242 is predicted to have a larger size than the normal Arg242, which may destabilize the protein structure. The hydrophobicity of the wild type and mutant residue differs. Hydrophobic interactions, either in the core of the protein or on the surface, will be lost. The homozygous c. 1015T>C (p. (Cys339Arg)) variant that was previously associated with autosomal recessive RP in a consanguineous family from Pakistan [53] is present three times in a homozygous state in ~60,000 'healthy' individuals in ExAC, which sheds doubt on its pathogenic character. Knockdown of znf513 in zebrafish model however suggested its role in retinal development and photoreceptor maintenance, as shown by the eye dysmorphology and photoreceptors loss in the knockdown model [53]. Possibly, the p. (Cys339Arg) variant shows reduced penetrance or is associated with variable age of onset of disease. 

Among the score 2 classified genes, we identified a homozygous truncating variant in TMED7. The role of TMED7 in protein trafficking and lack of an associated disease marks TMED7 as a putative target for further investigation. Similar to the other candidate genes, additional cases with variants in this gene would support its involvement in IRD. 

In conclusion, we have identified 11 candidate IRD-associated genes through data sharing. To prove that these candidate genes are true IRD causative genes, whole genome sequencing in these cases could further confirm their candidacy by exclusion of other genetic causes. Moreover, additional unrelated cases need to be identified and additional functional evidence can support causality. For this, we encourage a worldwide genotype-phenotype data sharing in order to boost the completion of the IRD genetic puzzle.
section_name: Discussion

section_content: The following are available online at www. mdpi. com/2073-4425/9/1/21/s1. Table S1 : Provided details on whole exome sequencing and linkage/homozygosity approach.
section_name: Supplementary Materials: