An evolutionary perspective into the role of kallikreins (KLKs) in male reproductive biology

PATRÍCIA ISABEL FERREIRA MARQUES AN EVOLUTIONARY PERSPECTIVE INTO THE ROLE OF KALLIKREINS (KLKs) IN MALE REPRODUCTIVE BIOLOGY Tese de Candidatura ao grau de Doutor em Ciências Biomédicas submetida ao Instituto de Ciências Biomédicas Abel Salazar da Universidade do Porto. Orientador – Doutora Susana Seixas Categoria – Investigadora Afiliação – Instituto de Investigação e Inovação em Saúde, Universidade do Porto (I3S); Instituto de Patologia e Imunologia Molecular da Universidade do Porto (Ipatimup). Coorientador – Doutor Victor Quesada Categoria – Investigador Afiliação – Departamento de Bioquímica y Biología Molecular de la Universidad de Oviedo Coorientador – Maria de Fátima Gärtner Categoria – Professora Catedrática Afiliação – Instituto de Ciências Biomédicas Abel Salazar da Universidade do Porto. iii Research work coordinated by: v Financiamento: Este trabalho foi financiado por Fundos FEDER através do Programa Operacional Factores de Competitividade – COMPETE e por Fundos Nacionais através da FCT – Fundação para a Ciência e a Tecnologia no âmbito do projeto FCOMP-01-0124-FEDER-028251 (Refª FCT: PTDC/BEX-GMG/0242/2012). Este trabalho foi ainda financiado pela FCT através da atribuição de uma bolsa individual de doutoramento (SFRH/BD/68940/2010). vii Ao abrigo do art.º 8º do Decreto-Lei n.º 388/70, fazem parte integrante desta dissertação os seguintes manuscritos já publicados, aceites para publicação ou em preparação: Marques PI, Bernardino R, Fernandes T. Nisc Comparative Sequencing Program, Green ED, Hurle B, Quesada V, Seixas S. 2012. Birth-and-Death of KLK3 and KLK2 in primates: evolution driven by reproductive biology. Genome Biol Evol. 4(12): 1331-8. Marques PI, Fonseca F, Sousa T, Santos P, Camilo V, Ferreira Z, Quesada V, Seixas S. 2015. Adaptive Evolution Favoring KLK4 Downregulation in East Asians. Mol Biol Evol. Epub ahead of print (DOI: 10.1093/molbev/msv199). Marques PI, Fonseca F, Carvalho AS, Puente DA, Damião I, Almeida V, Barros N, Barros A, Carvalho F, Mathiesen R, Quesada V, Seixas S. Rare and common variants in KLK and WFDC gene families and their implications into male infertility phenotypes. In preparation. Em cumprimento do disposto no referido Decreto-Lei, a candidata declara que participou na obtenção, análise e discussão dos resultados, bem como na elaboração das publicações, sob o nome Marques PI. “It is not our differences that divide us. It is our inability to recognize, accept, and celebrate those differences.” Audra Lorde xvi species the genomes are not yet fully assembled. On the left, a NCBI taxonomybased dendrogram shows the taxonomic classes and the evolutionary relationship among taxa. Data compiled from Pavlopoulou et al. 2010, Koumandou and Scorilas 2013 and Lundwall 2013. .......................................................................................17 Figure 4 –KLKs expression patterns in adult tissues. mRNA concentration for each KLK (as indicated in top row) and tissue. The color code at the bottom shows the levels of expression (from Shaw and Diamandis 2007). .........................................19 Figure 5 - Schematic representation of the semen liquefaction proteolytic cascade. (A) In normal physiologic conditions, KLKs are activated in the prostate through a zymogen activation cascade. KLK activation by other KLK is represented by straight arrows and auto-activation ability is illustrated by curved arrows. The pro-peptide is represented by the yellow rectangle. (B) Upon ejaculation, the sperm-rich epididymal fluid is mixed with prostatic fluids (including KLKs) along with secretions of the seminal vesicles (including SEMG1, SEMG2 and FN), forming the semen coagulum. SEMGs chelate Zn2+ ions, which leads to KLK reactivation and subsequent proteolysis of the SEMGs and FN, resulting in seminal coagulum liquefaction (adapted from Michael et al. 2006 and Prassas et al. 2015). ...............21 Figure 6 – Schematic representation of epidermis architecture and KLK proteolytic cascade in the skin. (A) The epidermis is organized in different layers mainly arranged by keratinocytes in different stages of differentiation. Keratinocytes are formed in the basal layer (stratum basale, SB) and begin to differentiate in the stratum spinosum (SS). This differentiation process occurs as keratinocytes migrate towards the skin surface. By the time these cells reach the stratum corneum (SC) they have already differentiated into corneocytes, cells filled with keratin and metabolically dead. (B) Pro-KLKs are secreted at the SG by lamellar granules (LG) of keratinocytes into SC interstices, where activation occurs by removal of the propeptide (yellow rectangle). Once active, KLKs cleave the corneodesmosome proteins, desmoglein 1 (DSG1), desmocollin 1 (DSC1) and corneodesmosin (CDSN), resulting in corneocyte shedding (skin desquamation). Several KLKs (KLK4, KLK5, KLK6 and KLK14) may also activate the protease-activated receptor2 (PAR-2), leading to inflammation, modulation of lipid-permeability barrier or melanosome transfer. The KLK activity in the skin is regulated by protease inhibitors, such as serine protease inhibitor Kazal-type 5 (SPINK5 or LEKTI), and by the epidermal pH gradient (adapted from Ovaere et al. 2009 and Prassas et al. 2015). ....................................................................................................................24 xvii Figure 7 – KLK4 in tooth enamel formation. The ameloblasts secrete a protein-rich matrix composed by amelogenin, enamelin and ameloblastin, as well as KLK4, MMP20 and DPP1 proteases. During the transitional and maturation stages, proKLK4 is secreted and activated by MMP20 and DPP1. Upon activation, KLK4 degrades the dental extracellular matrix proteins, allowing crystal growth in width and thickness, thus promoting enamel hardening (adapted from Prassas et al. 2015). ....................................................................................................................26 Chapter 3 Papers Paper I – Birth-and-Death of KLK3 and KLK2 in Primates: Evolution Driven by Reproductive Biology Figure 1 – Phylogenetic analysis of KLK2 and KLK3 in primates. (A) Phylogenetic tree showing primate divergence times (Hedges et al. 2006) and functional status of KLK2 and KLK3. The criteria to define a nonfunctional KLK gene were the identification of at least one disrupting mutation. Gray square indicates a duplication event. The ancestral KLK3 branch is indicated (ancKLK3). (B) Alignment of exons IV–V for KLK2 and KLK3 in Catarrhini. The corresponding human genomic positions for these regions are represented at the top. Positions conserved with Gorilla gorilla (left panel) or Nomascus leucogenys (right panel) are in orange. Nonconserved positions are in blue. Sites conserved in all species were omitted. . 38 Figure 2 – Positive selected sites in biologically relevant regions. (A) Human KLK2 three-dimensional model showing amino acid replacements predicted to be under positive selection (Q109, H177, and G210). (B) Human KLK3 three-dimensional model showing D207S substitution predicted to be under positive selection in the ancestral branch. The catalytic triad is represented in light blue (H65, D120, and S213) and the binding sites in orange (S228, G230, and D207 in KLK2 or S207 in KLK3). ................................................................................................................... 41 Figure 3 – Evolution of primate KLK2 andKLK3 related to mating factors. (A) Correlation of residual testis size (Anderson et al. 2004; Dixson and Anderson 2004; Wlasiuk and Nachman 2010) with the combined SEMG repeat units (JensenSeaman and Li 2003; Hurle et al. 2007). (B) Correlation between the number of xviii SEMG1 and SEMG2 repeat units (Jensen-Seaman and Li 2003; Hurle et al. 2007) and the presence of functional KLK2 and KLK3. *P<0.05. (C) Correlation between the mating system (Wlasiuk and Nachman 2010) and the presence of functional KLK2 and KLK3. UM, unimale; MM, multimale. *P<0.05. (●), monoandrous; (■), polyandrous; and (▲), ambiguous. ........................................................................ 42 Paper II – Adaptive Evolution Favoring KLK4 Downregulation in East Asians Figure 1 – Schematic representation of the human KLK gene cluster located at chromosome 19q13.3–13.4. Upper diagram shows the relative position of KLK genes. As depicted, the cluster includes 15 coding genes (black arrows) and one expressed pseudogene (gray arrow). The inset shows the KLK3–KLK5 region within the UCSC Genome Browser view for recombination maps from HapMap release 24 and UCSC gene transcripts. ................................................................. 48 Figure 2 – Sliding window of nucleotide diversity per base pair (x10-3) (A) and Tajima’s D (B) in the KLK3–KLK5 region in ASN (CHB+JPT), CEU, and YRI (solid, dashed, and dotted lines, respectively). Window size: 5,000 bp; increment: 1,000 bp. .............................................................................................. 51 Figure 3 – Genetic population differentiation (FST) analysis for KLK3–KLK5 locus of ASN versus non-ASN populations (A) and empirical rank FST scores based on global comparisons for CHB, CEU, and YRI (B). Genes location are delimited by open boxes. SNPs with significant FST P values (upper P < 0.05) or significant empirical rank scores (http://hsb.upf.edu/) are displayed in black. ......................... 52 Figure 4 – Ratio of intra-allelic diversity associated with the ancestral and derived alleles (iπA/iπD) plotted as a function of the DAF in the ASN (CHB+JPT) population. Black points: Candidate SNPs rs198968 and rs17800874. P < 0.05; solid line: 95% constant model; dashed line: 95% Laval mode (Laval et al. 2010); dotted line: 95% Gravel model (Gravel et al. 2011). ............................................... 53 Figure 5 – Signatures of natural selection at KLK3–KLK5 locus in human populations. (A) Worldwide estimated allele frequencies from 1000G data for variants rs1654556, rs198968, and rs17800874 in 14 human populations (Asia: CHB, JPT, and CHS—Southern Han Chinese in China; Africa: YRI and LWK— Luhya inWebuye, Kenya; Europe: CEU, GBR—British in England and Scotland, FIN—Finnish in Finland, IBS—Iberian populations in Spain, TSI—Toscani in Italia; Americas: ASW—African Ancestry in Southwest United States of America; CLM— xix Colombian in Medellin, Colombia, MXL—Mexican ancestry in Los Angeles, CA, PUR—Puerto Rican in Puerto Rico). (B) Schematic representation of ASN (CHB+JPT) haplotypes for KLK3–KLK5 region. Each line represents a haplotype and columns indicate polymorphic positions. Haplotypes are organized by different configurations of rs1654556, rs198968, and rs17800874 alleles. The relative positions of KLK genes are depicted by open arrows, the candidate SNPs by the filled arrows and the recombination hotspots (RH) are also shown. Ancestral alleles are represented in blue and derived alleles in orange. ........................................... 54 Figure 6 – In vitro validation of candidate variants by luciferase reporter assays. (A) pGL3 and pmirGLO constructs containing the ancestral (underlined) or the derived allele. The CNV alleles 105-bp deletion (Del105) and 67-bp insertion (Ins67) included in pmirGLO constructs are shown. Relative luciferase activity of variants rs198968 (B), rs17800874 (C), rs17800874+rs198968 (D), and rs1654556 (E) in LNCaP, HeLa, and AGS cell lines. Data are expressed as the mean ± standard error mean for at least three experiments.*P<0.05;**P<0.01; and ***P<0.0001. ..... 55 Paper III – Rare and common variants in KLK and WFDC gene families and their implications into semen hyperviscosity and other male infertility phenotypes Figure 1 – Minor allele frequencies (MAFs) from 1000 Genomes data vs. controls from pooled sequencing. Allele frequency estimates for 277 SNVs based on pooled sequencing from the control group were compared with the described European average frequencies from 1000 Genomes project phase III samples. r2 - correlation coefficient (r2 = 0.826). ......................................................................... 97 Figure 2 – Minor allele frequencies (MAFs) from pooled sequencing vs. Sanger sequencing. Estimated MAFs based on pooled sequencing is plotted against the actual frequencies as determined by individual Sanger sequencing for the surveyed regions. r2 - correlation coefficients (HV: r2 = 0.9725; NV: r2 = 0.9695; controls: r2 = 0.9509). The data from HV cases, NV cases and controls are represented in orange, blue and green, respectively. ................................................................... 97 Figure 3 – Structural characterization of the KLK low-frequency variants. (A) Alignment of the amino acid sequences of the variant kallikreins. Variant sites are framed in red. Complete conservation is shown in dark blue background, whereas partial conservation is shown on a light blue background. The catalytic serine is xx highlighted with a red arrow. (B) Mapping of variant sites on a kallikrein structure. The overall structure is depicted as a green ribbon. Variant sites are shown as sticks. The catalytic triad and the second SS6 cysteine are shown as lines. ......... 98 Figure 4 – Relative abundance of KLK3 p.S210W variant in seminal plasma. Spectral counts for p.S210W residue in two heterozygous (Het_1, Het_2) and of 33 homozygous (Hz) individuals. Total spectral counts are shown for Het_1 and Het_2 individuals, and the mean of spectral counts are displayed for Het_1+2 and Hz. ... 99 Chapter 4 Final Discussion Figure 8 – Worldwide estimated haplotype frequencies defined by rs1654556, rs198968 and rs17800874 according to 1000G phase III data for African, European, South Asia, East Asia and American populations. For each continental region the most common haplotypes are shown. In Africa the ancestral haplotype is also displayed. Ancestral and derived alleles are represented in blue and orange, respectively. ..................................................................................... 107 Appendices Appendix A – Supplementary Material Paper I Figure S1 - KLK3-KLK2 gene fusion event in Gorilla gorilla, Nomascus leucogenys and Hylobates sp. Schematic representation of G. gorilla (A) or N. leucogenys (B) genomic sequence alignments with the Homo sapiens reference sequence (KLK3 to KLK4). BlastN hits are represented as boxes joined with a line. Lighter lines indicate a non-optimal hit in one of the regions. Insertions and deletions cause a lack of correspondence between sequences. (C) Gene fusion event confirmed in G. gorilla and Hylobates sp. by PCR assay with genespecific primers for KLK3 (exon 4) and for KLK2 (exon 5). The gene fusion product was confirmed in both species by sequencing of the resulting amplicons. ............................................................ 149 Figure S2 - Genomic sequence alignments of the orthologous genomic fragment spanning KLK1to KLK4 in Colobus guereza and Homo sapiens (green and red, respectively). BlastN hits are represented as boxes joined with a line. Lighter lines xxi indicate a non-optimal hit in one of the species. Insertions and deletions cause lack of correspondence between sequences. .............................................................. 150 Figure S3 – KLK2 protein alignment identifying deleterious mutations. ( ) Start codon; ( ) Catalytic triad residues; ( ) Activation site; ( ) Frameshift; ( ) Premature STOP codons; (?) Missing data. .......................................................................... 151 Figure S4 – Evolution of primate KLK2 and KLK3 related to reproductive traits. A) Correlation between the presence of functional KLK2 and KLK3 with semen coagulation rating. Semen coagulation is rated on a four-point scale (Dixson and Anderson 2002), with 1 reflecting no coagulation and 4 reflecting the production of a solid copulatory plug. B) Correlation of residual testis size (Anderson et al 2004; Dixson and Anderson 2004; Wlasiuk and Nachman 2010) with the presence of functional KLK2 and KLK3. .................................................................................. 152 Appendix B – Supplementary Material Paper II Figure S1 – Selection statistics for KLK3-KLK5 locus. (A) Cross-Population Extended Haplotype Homozygosity (XP-EHH) plot from HGDP data for different continental populations as indicated by different color lines (http://hgdp.uchicago.edu/cgibin/gbrowse/HGDP/). East Asia is represented in green, South Asia in black, Europe in orange, Mideast in blue, Oceania in turquoise, America in yellow, Bantu in red and non-Bantu African populations in pink and purple. (B) 1000 Genomes Selection Browser view. Statistic tracks for pairwise FST for CHB vs. CEU, YRI vs. CHB and CEU vs. YRI, FST Global (CHB, CEU and YRI), integrated haplotype score (iHS) for CHB, cross-population extended haplotype homozygosity (XP-EHH) for CHB vs. CEU and YRI vs. CHB, and cross-population composite likelihood ratio (XP-CLR) for CHB vs. CEU and YRI vs. CHB. The statistics are presented as – log10 of empirical ranked scores (http://hsb.upf.edu/). ......................................... 177 Figure S2 – Genetic population differentiation (FST) analysis for KLK3-KLK5 locus of ASN vs. CEU, ASN vs. YRI and CEU vs. YRI populations. Genes’ location is delimited by open boxes. SNPs with significant FST P-values (upper P < 0.05) are displayed in blue, green and red for ASN vs. CEU, ASN vs. YRI and CEU vs. YRI comparisons, respectively. ................................................................................... 178 Figure S3 – Linkage disequilibrium plot of 1000G phase I data for KLK3-KLK5 region in Asians. The image was generated using Haploview 4.2 software. The triangular units represent haplotype blocks as defined by Gabriel et al. 2002. The xxii degree of LD between pair of markers is indicated by the |D’| statistic (|D’| = 1, bright red; |D’| < 1, shades of red). The relative positions of KLK genes are depicted by open arrows, and the relative positions of the recombination hotspots are also shown. ................................................................................................................. 178 Figure S4 - Schematic representation of KLK3-KLK5 landscape using UCSC Genome Browser. Reference genes, DNase hypersensitivity and chromatin state segmentation from ENCODE are shown in the upper image. The insets display in detail the KLK4 locus and the putative enhancer within the intergenic region between KLK4 and KLK5. The SNPs rs1654556, rs198968 and rs17800874 are highlighted by red circles. ................................................................................... 179 Figure S5 – Worldwide allele frequencies from HGDP data for rs198968 and 17800874 SNPs as inferred by fastPHASE (adapted from http://hgdp.uchicago.edu/cgi-bin/gbrowse/HGDP/). (A) Frequencies of rs198968 located in intron I of KLK4. (B) Frequencies of rs17800874 located in a putative enhancer in the intergenic region between KLK4 and KLK5. ................................ 180 Figure S6 – Extended haplotype homozygosity (EHH) statistic for ASN (CHB+JPT) sample using 1000G data. Plots of EHH over genetic distance for the largest nonoverlapping cores encompassing rs1654556 (A), rs198968 (B) or rs17800874 (C) variants. Core haplotype sequences are indicated below EHH plots and candidate variants underlined. ............................................................................................. 181 Figure S7 – Plots of KLK4 expression for rs198968 (A) and rs17800874 (B) quantitative trait loci (eQTL) in prostate tissues from GTEx data (http://www.gtexportal.org/home/). The corresponding genotypes are indicated in parenthesis and the number of samples and P-values are shown........................ 181 Figure S8 – Tissue expression of KLK2, KLK3, KLK4, KLK5 genes and KLKP1 pseudogene. Multiplex PCRs carried out in a cDNA panel from human healthy organs, each one including a minimum of three donor’s pool. GAPDH or SERPINA1 fragments were used as internal controls. ............................................................ 182 Appendix C – Supplementary Material Paper III Figure S1 – Flow-chart of the strategy used to detect rare and common variants on KLK and WFDC clusters associated with male infertility. Using a DNA pooled sample approach and a high-throughput sequencing strategy, we detected in phase I 456 SNVs based on stringent filtering criteria. We then performed genotyping xxiii validation of 3 SNVs and 7 gene regions in phase II, using the same samples as in phase I. In phase III, we extended the analysis of the most promising SNVs to a further 138 controls and 95 infertility cases to allow a combined analysis of 217 controls and 238 cases. ....................................................................................... 199 Figure S2 – Schematic representation of the human KLK and WFDC gene clusters using UCSC Genome Browser. (A) The human KLK cluster is located on chromosome 19q13.3-13.4 and includes 15 coding genes and one expressed pseudogene. (B) The human WFDC cluster located is on chromosome 20q13 and its genes are organized into two subloci (centromeric and telomeric, WFDC-CEN and WFDC-TEL, respectively), separated by 215 kb of unrelated sequence. Amplicons generated for the pilot survey, reference genes, H3K4Me1 Mark, DNase hypersensitivity and transcription factor CHIP-seq from ENCODE are shown. ..... 200 Figure S3 – WFDCs minor allele frequencies (MAFs) from 1000 Genomes data and control pooled sequencing in repetitive regions. Allele frequency estimates obtained in pooled sequencing for the control group (black) and the described frequencies from the combined European populations from 1000 Genomes project phase III (orange). ............................................................................................... 201 Figure S4 – KLKs minor allele frequencies (MAFs) from 1000 Genomes data and control pooled sequencing in repetitive regions. Allele frequency estimates obtained in pooled sequencing for the control group (black) and the described frequencies from the combined European populations from 1000 Genomes project phase III (orange). ............................................................................................... 202 Figure S5 – Schematic representation of the human KLK7 landscape using UCSC Genome Browser. Amplicons generated for the pilot survey, reference genes, H3K4Me1 Mark, DNase hypersensitivity, transcription factor CHIP-seq and chromatic state segmentation from ENCODE are shown in the upper image. The inset displays in detail the intron V of KLK7 in which rs1654526 is located (highlighted by red circle). .................................................................................... 203 Figure S6 – Alignment of the kallikrein protein sequences. Complete conservation is shown in dark blue background, whereas partial conservation is shown on a light blue background. The catalytic residues are framed in red. Variant sites are indicated by arrows. The equivalent variants 131_KLK3 and 138_KLK14 are highlighted in pink. ............................................................................................... 204 xxv Tables List Chapter 1 General Introduction Table 1 - Common approaches used to detect selection (adapted from Vitti et al. 2013). ....................................................................................................................11 Table 2 – Semen quality nomenclature according to WHO 1999. ...............................22 Chapter 3 Papers Paper I – Birth-and-Death of KLK3 and KLK2 in Primates: Evolution Driven by Reproductive Biology Table 1 – Identified KLK2 Deleterious Mutations. ........................................................ 39 Table 2 – Parameter Estimates and Likelihood Scores under Different Branch Models. ................................................................................................................. 40 Table 3 – Model Comparisons of Variable o Ratios among Sites. .............................. 41 Paper II – Adaptive Evolution Favoring KLK4 Downregulation in East Asians Table 1 – Summary Statistics of KLK3–KLK5 Population Variation from 1000G Data…. .................................................................................................................. 50 Paper III – Rare and common variants in KLK and WFDC gene families and their implications into semen hyperviscosity and other male infertility phenotypes Table 1 – Burden tests for KLK and WFDC low-frequency variants. .......................... 93 Table 2 – Low-frequency variants surveyed in phase II. ............................................. 94 xxxii KLK11 - Kallikrein11 (also known as hippostasin or serine protease 20) KLK12 - Kallikrein 12 (also known as KLK-L5) KLK13 - Kallikrein 13 (also known as KLK-L4) KLK14 - Kallikrein 14 (also known as KLK-L6) KLK15 - Kallikrein 15 (also known as prostinogen) KLK1E2 - Equus caballus glandular kallikrein precursor KLK1P - Canis lupus familiaris kallikrein 1 pseudogene KLK2 - Kallikrein 2 (also known as human glandular kallikrein-1 or tissue kallikrein-2) KLK3 - Kallikrein 3 (also known as prostate-specific antigen) KLK4 - Kallikrein 4 (also known as KLK-L1, enamel matrix serine protease 1 or prostase serine protease 17) KLK5 - Kallikrein 5 (also known as KLK-L2 or stratum corneum tryptic enzyme) KLK6 - Kallikrein 6 (also known as neurosine) KLK7 - Kallikrein 7 (also known as stratum corneum chymotryptic enzyme) KLK8 - Kallikrein 8 (also known as neuropsin) KLK9 - Kallikrein 9 (also known as KLK-L3) KLKP1 - Kallikrein pseudogene 1 KRT77 - Keratin 77, type II L L LCT - lactase LD - Linkage disequilibrium LEKTI -Llymphoepithelial kazal type inhibitor LG - Lamellar granules LRH - Long-range haplotype LWK - Luhya inWebuye in Kenya xxxiii M M MAF - Minor allele frequency ME2 - Malic enzyme 2, NAD(+)-dependent, mitochondrial ME3 - Malic enzyme 3, NADP(+)-dependent, mitochondrial miRNA - micro ribonucleic acid MK - McDonald-Kreitman mL - Mililiter MM - Multimale MMP20 - Matrix metalloproteinase-20 mRNA - messenger ribonucleic acid MS – mass spectrometry MTRR - 5-methyltetrahydrofolate-homocysteine methyltransferase reductase MXL - Mexican ancestry in Los Angeles, CA mya - Million years ago N N N/A – Not applicable NGS – Next-generation sequencing NS - Netherton syndrome NV – Non-hyperviscosity O O OCA2 - Oculocutaneous albinism II P P PAR-2 - Protease-activated receptor-2 PCR - Polymerase chain reaction PI3 - Peptidase inhibitor 3, skin-derived (also known as WFDC14 or ELAFIN) xxxiv PLRP2 - Pancreatic lipase-related protein 2 PSA - Prostate-specific antigen (also known as KLK3) PUR - Puerto Rican in Puerto Rico R R r2 – linear correlation coefficient REHH - Relative extended haplotype homozygosity RH - Recombination hotspot RNA - Ribonucleic acid RPTOR - Regulatory associated protein of MTOR, complex 1 RT-PCR - Reverse transcription polymerase chain reaction S S SB - Stratum basale SC - Stratum corneum SDS-PAGE - Sodium dodecyl sulfate poly-acrylamide gel electrophoresis SEMG - Semenogelins SEMG1 - Semenogelin 1 SEMG2 - Semenogelin 2 SERPINA5 - Serpin peptidase inhibitor, clade A (alpha-1 antiproteinase, antitrypsin), member 5 SERPINB11 - Serpin peptidase inhibitor, clade B (ovalbumin), member 11 SERPINB3 - Serpin peptidase inhibitor, clade B (ovalbumin), member 3 SERPINB4 - Serpin peptidase inhibitor, clade B (ovalbumin), member 4 SFS - Site frequency spectrum SG - Stratum corneum SIGLEC5 - Sialic acid binding Ig-like lectin 5 SIGLEC6 - Sialic acid binding Ig-like lectin 6 xxxv SLC24A5 - Solute carrier family 24 (sodium/potassium/calcium exchanger), member 5 SLC45A2 - Solute carrier family 45, member 2 SLPI - Secretory leukocyte peptidase inhibitor (also known as WFDC4) SMIPS - Somatic Mutation Identification in Pooled Samples SNP - Single nucleotide polymorphism SNV - Single nucleotide variant SPINK5 - Serine protease inhibitor Kazal-type 5 SPINLW1 - Serine protease inhibitor-like with Kunitz and WAP domains 1 (also known as EPPIN) SRY - Sex determining region Y SS - stratum spinosum T T TGM4 - Transglutaminase 4 TRH - thyrotropin-releasing hormone TYRP1 – Tyrosinase-related protein 1 TSI - Toscani in Italia U U UM - Unimale UTR - Untranslated region UV-light - Ultraviolet light V V VEP - Variant Effect Predictor W W WFDC - whey acidic protein four-disulfide core domain WHO - World Health Organization xxxvi WNT10A - Wingless-type MMTV integration site family, member 10A X X XP-CLR - Cross-population composite likelihood ratio XP-EHH - Cross-population extended haplotype homozygosity Y Y YRI - Yoruba in Ibadan, Nigeria Z Z Zn2+ - Zinc ZP - zona pellucida ZP2 - Zona pellucida glycoprotein 2 ZP3 - Zona pellucida glycoprotein 3 1 Abstract The advances in sequencing technologies have greatly contributed to the increasing number of available genomes for different species, as well as to the development of detailed catalogs of human genetic variation. Likewise, these represent important tools for detecting footprints of natural selection acting on different timescales and help providing a better understanding of the molecular basis of current patterns of human disease susceptibility. Among the myriad of targets of natural selection, genes involved in reproductive functions are particularly relevant since they are at the frontline of the individual fitness. The kallikrein (KLK) gene family encodes 15 serine proteases, often co-expressed in a wide variety of tissues and in many biological fluids, including in the seminal plasma. In this sense, KLKs are known to modulate key physiological processes through complex proteolytic cascades, in which family members act at different hierarchical levels. In the particular case of the semen liquefaction cascade, KLKs are involved in the hydrolysis of the major seminal structural proteins, the semenogelins (SEMGs), resulting in sperm release in the vaginal cavity. Previous studies provided evidence of KLKs substrates being preferred targets of natural selection through mechanisms linked to male fertility and sperm competition. Also, the loss of KLK2 in some primate species was found to correlate to different semen physiologies. In addition, aberrant expression of most human KLK members was reported in individuals with abnormal semen parameters. The three independent studies presented in the present work were all based on the central hypothesis that KLK genes might have been targeted by natural selection and that their genetic variation may underlie both beneficial and disease phenotypes. The analysis of the selective pressures acting on the KLK cluster was performed in two stages. The first was based on comparative and phylogenetic approaches, centered in KLK2 and KLK3, for a total of 22 primate species. The second enclosed a comprehensive evaluation of the 1000 Genomes phase I data, in order to characterize a potential signature of natural selection among KLK genes in Asian populations. At the interspecific level, this study supported the origin of KLK3 in Catarrhini through an event of KLK2 duplication and functional divergence of KLK3 towards a different substrate specificity and it unraveled an intricate evolutionary dynamics of KLK2 and KLK3 correlated to semenogelin gene structure, primate mating system and semen coagulation rates. On the other hand, in human populations, a complex signature of recent positive selection in East Asians was disclosed and characterized by a high frequency haplotype defined by three variants (rs1654456_G, rs198968_T and rs17800874_A) acting synergistically to promote KLK4 2 downregulation, which may be connected to tooth and epidermal features typical of these populations rather than with reproductive functions. In the last part of this work, the variation of KLK genes, their substrates (semenogelins) and potential inhibitors of the whey acidic protein four-disulfide core domain (WFDC) locus was assessed in the scope of male infertility. This analysis revealed a higher burden of functional low-frequency variants in cases than in controls, independently of the considered infertility phenotype, but solely for the KLK cluster. Furthermore, it was possible to identify a significantly increased risk of semen hyperviscosity and asthenozoospermia associated with the variants rs61742847 (KLK12, p.P34L) and rs147894843 (SEMG1, p.G400D), respectively. In addition, other 12 nucleotide variants were also overrepresented in cases but, due to the still limited number of samples used in the genotype surveying, those were not statistically significant. Conversely, a decreased risk of hyperviscosity and oligozoospermia was observed for rs1654526 in KLK7 and for a copy number variation in SEMG1, correspondingly. Altogether, this work further supports that genes involved in proteolysis and reproductive biology, such as KLK genes, were targets of natural selection in the short and the larger timescales of human and primate evolution, respectively, and that KLK genetic variation may also underlie deleterious variants possible contributing to a lower reproductive fitness in humans. 3 Resumo Os progressos nas tecnologias de sequenciação têm sido um forte contributo para o aumento do número de genomas disponíveis para diferentes espécies, bem como para o desenvolvimento de catálogos detalhados da variação genética humana. Deste modo, estes representam importantes ferramentas para a deteção de evidências de seleção natural em diferentes escalas temporais e ajudam ainda a obter um maior conhecimento das bases moleculares dos atuais padrões de suscetibilidade a doenças humanas. Entre os vastos alvos de seleção natural, os genes envolvidos em funções reprodutivas são particularmente relevantes, uma vez que estão na vanguarda da aptidão de um indivíduo. A família de genes das calicreínas (KLK) codifica 15 proteases de serina, frequentemente co-expressas numa grande variedade de tecidos e em muitos fluidos biológicos, incluindo no plasma seminal. Neste sentido, as KLKs são conhecidas por modular importantes processos fisiológicos através de cascatas proteolíticas complexas, em que membros desta família atuam em diferentes níveis hierárquicos. No caso particular da cascata de liquefação do sémen, as KLKs estão envolvidas na hidrólise das principais proteínas estruturais, as semenogelinas (SEMGs), resultando na libertação dos espermatozóides dentro da cavidade vaginal. Estudos anteriores mostraram que os substratos das KLKs foram alvos de seleção natural através de mecanismos relacionados com a fertilidade masculina e competição espérmica. Também a perda da KLK2 em algumas espécies de primatas está correlacionada com diferentes fisiologias do sémen. Além disso, foi reportada uma expressão aberrante de KLKs humanas em indivíduos com parâmetros seminais anormais. Os três estudos independentes apresentados no presente trabalho foram baseados na hipótese central de que genes das KLKs podem ter sido alvo de seleção natural e que a sua variação genética pode contribuir para fenótipos benéficos e de doença. A análise das pressões seletivas no agrupamento das KLKs foi realizada em duas etapas. A primeira baseou-se em abordagens comparativas e filogenéticas, centradas na KLK2 e na KLK3, para um total de 22 espécies de primatas. A segunda compreendeu uma avaliação abrangente dos dados da fase I dos 1000 genomas, com a finalidade de caracterizar uma potencial assinatura de seleção natural entre genes das KLKs em populações asiáticas. A nível interespecífico, este estudo suporta a origem da KLK3 nos Catarríneos através de um evento de duplicação da KLK2 e divergência funcional da KLK3 para uma diferente especificidade de substrato e desvenda uma intrincada dinâmica evolutiva da KLK2 e da KLK3 correlacionada com a estrutura génica das semenogelinas, sistemas de acasalamento nos primatas e rácios de coagulação do 4 sémen. Por outro lado, nas populações humanas, foi revelada uma complexa assinatura de seleção positiva recente nos asiáticos de leste, caracterizada por um haplótipo de elevada frequência, definido por três variantes (rs1654456_G, rs198968_T e rs17800874_A) que atuam sinergicamente para promover uma redução nos níveis de KLK4, o que pode estar relacionado com características dentárias e epidérmicas típicas dessas populações e não com funções reprodutivas. Na última parte deste trabalho, foi avaliada no âmbito da infertilidade masculina a variação dos genes das KLKs, dos seus substratos (semenogelinas) e potenciais inibidores do locus whey acidic protein four-disulfide core domain (WFDC). Esta análise revelou uma maior carga de variantes funcionais de baixa frequência nos casos do que nos controlos, exclusivamente no agrupamento das KLK e independentemente do fenótipo de infertilidade considerado. Além disso, foi possível identificar um maior risco significativo de hiperviscosidade do sémen e de astenozoospermia associado com os variantes rs61742847 (KLK12, p.P34L) e rs147894843 (SEMG1, p.G400D), respetivamente. Adicionalmente, outros 12 variantes nucleotídicos encontravam-se também sobre-representados em casos mas, devido ao limitado número de amostras utilizadas na genótipagem, estes não atingiram significância estatística. Por outro lado, foi observada uma diminuição do risco de hiperviscosidade e de oligozoospermia para o variante rs1654526 na KLK7 e para uma variação do número de cópias (CNV) na SEMG1, correspondentemente. De um modo geral, este trabalho suporta que genes envolvidos na proteólise e na biologia reprodutiva, tais como os genes das KLKs, foram alvos de seleção natural ao longo da escala evolutiva dos humanos e dos primatas. Este estudo suporta ainda que a variação genética das KLKs pode também conter variantes deletérios que poderão contribuir para um menor fitness reprodutivo nos humanos. Chapter 1 General Introduction Chapter 1 | General Introduction 12 Nielsen 2002; Bielawski and Yang 2003; Yang et al. 2005). These tools have been used in a wide variety of cases, for example, to reveal the presence of positive selection on genes involved in language and speech (FOXP2), reproduction (ZP2, ZP3 and ADAM2) and immunity (SIGLEC5, SIGLEC6, SERPINB3 and SERPINB4) (Swanson et al. 2001; Enard et al. 2002; Zhang et al. 2002; Kosiol et al. 2008; Gomes et al. 2014). Recently, there is a growing body of evidence that regions of the human genome targeted by positive selection may also be associated with human disease. One possible explanation for this correlation is a change in the selective forces acting on human populations; alleles which were once beneficial may now be deleterious because the present environmental conditions are not the same as in the past (“thrifty genotype” hypothesis) (Neel 1962; Neel et al. 1998). This principle is well illustrated by variants predisposing to type-II diabetes. Several generations ago, fat accumulation was advantageous because it would largely raise the chances of survival in times of famine and, therefore, variants associated with this trait were positively selected. However, in modern societies, the current settled lifestyle and the ample food availability no longer render this characteristic beneficial and variants that were once valuable become maladaptive by increasing the risk of their carrier developing type-II diabetes (Di Rienzo and Hudson 2005; Di Rienzo 2006; Helgason et al. 2007; Crespi 2010; Vasseur and Quintana-Murci 2013). On the other hand, selective events could end in an antagonist pleiotropy, a phenomenon in which selection results in adaptation in one trait, or early in life, and in deleterious effects in other contexts or later in the lifespan (Clark and Swanson 2005; Nielsen et al. 2005; Corbo et al. 2008; Crespi 2010; Vasseur and Quintana-Murci 2013). Furthermore, when an advantageous variant rises in frequency during the selective sweep, the hitchhiking effect can drive disease-causing alleles to high frequency as well (Shiina et al. 2006; Huff et al. 2012). Therefore, there may be a close link between selection and disease. Thus, the combination of comparative and population genetic tools with functional information can be useful to exploit such disease candidate loci. Overall, comparative genomics, human genome-wide scans (GWS) and candidate gene approaches have identified many potential examples of selection targets and verified that certain gene ontology categories, including sensory perception, dietary changes, immunity and host-pathogen interactions, reproduction and proteolysis, were enriched in genes under selection (Bustamante et al. 2005; Sabeti et al. 2006; Kosiol et al. 2008; Akey 2009; Bustamante and Ramachandran 2009; Grossman et al. 2013). In this context, the kallikrein (KLK) gene family represents a remarkable case for the study of evolution of proteolytic genes with implications in different biological processes including in reproductive biology and in human health and disease. Chapter 1 | General Introduction 13 2. The Kallikrein (KLK) locus The first kallikrein was identified in the 1930s as an abundant protease in the pancreas and, therefore, was named tissue kallikrein, based on the Greek word for pancreas “kallikreas”. Later work from independent groups led to the identification and characterization of 14 additional genes that now comprise the KLK gene family (Lilja 1985; Riegman et al. 1992; Gan et al. 2000; Yousef et al. 2000; Clements et al. 2001). As the members KLK1, KLK2 and KLK3 were the first ones described, they became known as “classical kallikreins” and KLK1 remained as the prototypical kallikrein gene. In order to avoid confusion and harmonize the terminology of these genes, an official nomenclature system was proposed by the Human Genome Organization (HUGO). According to this nomenclature, KLK1 is called kallikrein-1, whereas every other member of the family is termed kallikrein-related peptidase (Lundwall et al. 2006a). However, all these molecules are often simply named as kallikreins. Since their discovery, it became evident that KLKs are involved in many physiological processes and in several pathophysiologic conditions, as described below. 2.1. Structure and organization The human KLK cluster, located at chromosome 19q13.3-13.4, spans over 265 kb and includes 15 paralogue genes coding for trypsinor chymotrypsin-like serine proteases (KLK1 to KLK15) and a transcribed pseudogene (KLKP1) (Yousef et al. 2004; Lu et al. 2006; Kaushal et al. 2008; Lundwall and Brattsand 2008; Prassas et al. 2015). The intergenic spacing between KLK genes is variable and ranges from approximately 1.5 to 32.5 kb, in which the smaller intergenic region is located between KLK1 and KLK15 and the largest one between KLK4 and KLK5. Typically, genes extend from 4.4 to 10.5 kb, depending mainly on intron sizing, and, except for KLK2 and KLK3, all of them are transcribed in the reverse chromosome strand. Moreover, all KLK genes display a common organization in five coding exons with variable 5’ and 3’ untranslated region (UTR) structures (Figure 2A) (Obiezu and Diamandis 2005; Lawrence et al. 2010). Furthermore, all KLK genes have at least one additional transcript resulting from alternative splicing, which in most instances are non-coding RNAs or result in shorter proteins with no protease activity (Kurlender et al. 2005; Koumandou and Scorilas 2013). Nonetheless, the function of these transcripts remains poorly understood and it has been suggested that they may mediate the signaling either at mRNA or at protein level (Koumandou and Scorilas 2013). The only exception to these structural features is KLKP1 Chapter 1 | General Introduction 14 which, in spite of having five transcribed exons, the in silico translation of the four mRNA transcripts indicate that none of them encode a serine protease and the two largest putative proteins are generated from a single exon (Yousef et al. 2004; Lu et al. 2006; Kaushal et al. 2008). Figure 2 – Genomic and proteomic structure of KLK proteases. (A) The human KLK gene cluster is located at chromosome 19q13.3-13.4. Arrows show the relative position and the transcription orientation for the 15 coding genes and expressed pseudogene. In the mRNA scheme, the boxes and lines represent exons and introns, respectively. KLK proteins are expressed as pre-pro-enzymes, in which the pre-domain is required for intracellular trafficking and the pro-domain must be cleaved in order to generate a mature KLK. (B) The structure of a mature KLK based on the crystal structure of KLK1. The catalytic triad residues are shown in green. The position of the kallikrein loop is also shown. The amino acids that are identical among all kallikreins are in red, whereas those that are conserved in at least eight kallikreins are in purple. Nonconserved amino acids are in blue (adapted from Lawrence, Lai, and Clements 2010 and Prassas et al. 2015). Chapter 1 | General Introduction 15 Kallikreins are encoded as single-chain pre-pro-enzymes (Figure 2A) of 244 to 293 residues long, sharing about 40-80% of protein identity, as well as a conserved catalytic triad of histidine (H57), aspartic acid (D102) and serine (S195) residues (standard chymotrypsin numbering) (Lundwall and Brattsand 2008). Besides the catalytic triad, amino acids involved in protein folding are also highly conserved among KLKs, whereas residues associated with substrate specificity are generally more divergent (Figure 2B). Once synthesized KLKs are directed to the endoplasmic reticulum by the signal peptide (the pre-domain), which is composed by 16 to 33 N-terminal amino acids, this pre-peptide is later cleaved in the secretory pathway yielding an enzymatically inactive pro-KLK. These molecules (zymogens) only become active upon the proteolytical removal of the pro-domain, which allows a conformational rearrangement of the protein threedimensional structure by opening the KLK catalytic fissure. In most cases, the KLK prodomain is cleaved after an arginine or lysine residue, indicating that they are activated by proteases with trypsin-like specificity including other KLKs, or in some instances by themselves (auto-activation) (Vaisanen et al. 1999; Brattsand et al. 2005; Michael et al. 2005; Memari et al. 2007; Yoon et al. 2007; Yoon et al. 2009; Lawrence et al. 2010). The exception to this rule is KLK4, which is cleaved after a glutamine residue and instead is activated by matrix metalloproteinase-20 (MMP20) and dipeptidyl peptidase 1 (DPP1) (Ryu et al. 2002; Tye et al. 2009; Yamakoshi et al. 2013). So far, six human KLK structures have been solved (KLK1 and KLK3 to KLK7) and, according to those crystallographic models, KLKs are folded into two hydrophobic interacting sub-domains, each one comprising a six-stranded β-barrel and a α-helix, with the catalytic triad located at the interface between the two sub-domains (Figure 2B) (Bernett et al. 2002; Gomis-Ruth et al. 2002; Laxmikanthan et al. 2005; Debela et al. 2006; Debela et al. 2007a; Debela et al. 2007b; Debela et al. 2008; Menez et al. 2008). KLK substrate specificity depends on the residue that lies at the base of the substrate binding pocket, and it is further refined by the composition of the eight loops that surround the active site (Debela et al. 2008; Lawrence et al. 2010). In addition, KLK1, KLK2 and KLK3 structures also contain an 11-amino acid insertion known as the “kallikrein loop”, which is believed to play an important role in the substrate and inhibitor specificity (Yousef and Diamandis 2001; Borgono and Diamandis 2004; Laxmikanthan et al. 2005; Lundwall and Brattsand 2008). Chapter 1 | General Introduction 16 2.2. Phylogenetic evolution The structure of KLK genes and their organization within a single syntenic locus suggest that all members of this family have evolved from a common ancestor by a series of gene duplication events that occurred at different moments of vertebrate evolution (Figure 3). First, it was suggested that the origin of the KLK family arose before the marsupial-placental split, approximately 125-175 million years ago (mya) (Elliott et al. 2006). However, the growing number of genomes available in public databases have enabled a better resolution of the KLK cluster evolutionary history, tracing its origins to a tetrapod lineage approximately 330 mya (Pavlopoulou et al. 2010; Koumandou and Scorilas 2013; Kawasaki et al. 2014). Moreover, the identification of 11 marsupial orthologs (KLK5 to KLK15) places the majority of duplication events prior to the separation of marsupial (Metatheria) and placental mammals (Eutheria) (Pavlopoulou et al. 2010; Koumandou and Scorilas 2013; Kawasaki et al. 2014). Consistently, several KLK orthologs were described for the platypus genome (Ornithorhynchus anatinus), but their tandem organization in a single syntenic cluster could not be assessed due to the lack of contiguous genomic segments and a still incomplete genome assembly (Pavlopoulou et al. 2010; Koumandou and Scorilas 2013; Lundwall 2013) (Figure 3). Although the phylogenetic relationships among KLK5 to KLK15 genes are not yet well resolved, most studies done so far seem to agree that a single duplication encompassing KLK9 and KLK10 yielded KLK11 and KLK12, or vice-versa (Olsson et al. 2004; Elliott et al. 2006; Lundwall et al. 2006b; Pavlopoulou et al. 2010; Lundwall 2013). Conversely, the events that generated KLK2, KLK3 and KLK4 are thought to have happened only after the placental mammals split (Elliott et al. 2006; Pavlopoulou et al. 2010; Koumandou and Scorilas 2013; Lundwall 2013). Initially, the absence of KLK4 in elephant, hyrax, tenrec and armadillo genomes, together with phylogenetic data, suggested that this gene had arisen by a duplication of KLK5 in the Boreoeutheria lineage (Elliott et al. 2006; Pavlopoulou et al. 2010; Koumandou and Scorilas 2013; Lundwall 2013). Still, the recent identification of a fragmented KLK4 pseudogene in the genome of cape golden mole (Chrysochloris asiatica; Afrotheria) points to an earlier duplication event that most likely occurred in a common ancestor of placental mammals (Kawasaki et al. 2014). Further comparative analyses show even more complex histories for the classical KLKs, which display several events of gene gain and loss. For instance, both the dog (Canis lupus familiaris) and horse (Equus caballus) genomes present an additional KLK1 homolog but, while in the horse this gene remained functional (KLK1E2), in the dog it degenerated into a pseudogene (KLK1P). On the other hand, in rodents, a remarkable Chapter 1 | General Introduction 17 Figure 3 – Schematic representation of KLK genes in different species. The arrows specify the direction of transcription and known pseudogenes are indicated in red. Loci are not drawn to scale and bars do not represent chromosomes, as for many species the genomes are not yet fully assembled. On the left, a NCBI taxonomy-based dendrogram shows the taxonomic classes and the evolutionary relationship among taxa. Data compiled from Pavlopoulou et al. 2010, Koumandou and Scorilas 2013 and Lundwall 2013. Chapter 1 | General Introduction 18 number of gene duplications generated an expanded locus with 13 and 10 Klk1 paralogs in mouse (Mus musculus) and rat (Rattus norvegicus), respectively, as well as several pseudogenes in mouse (Puente et al. 2003; Puente and López-Otín 2004; Pavlopoulou et al. 2010; Koumandou and Scorilas 2013; Lundwall 2013). Moreover, a KLK1 duplication leading to KLK2 seems to have occurred early in the mammal tree (Boreoeutheria), however, this appears to have been deleted or silenced into a pseudogene in many species. Finally, in the primate lineage, in a common ancestor of Old Word Monkeys (Catarrhini), a recent duplication of KLK2 generated KLK3 (Elliott et al. 2006; Pavlopoulou et al. 2010) (Figure 3). 2.3. Biological importance in human health and disease Kallikreins are often co-expressed in a wide variety of tissues (Figure 4) and found in many biological fluids. Accordingly, KLKs have been implicated in a broad range of physiological functions and proteolytic cascades, including semen liquefaction, skin desquamation, tooth enamel formation, neural plasticity and regulation of blood pressure (Klokk et al. 2006; Shaw and Diamandis 2007; Lundwall and Brattsand 2008; Lawrence et al. 2010; Prassas et al. 2015). Moreover, the activity of KLKs is generally regulated by fine-tuned processes and spatial-temporal modifications of KLK activity have been associated with pathological conditions such as psoriasis, atopic dermatitis, hypertension, diabetes, neurodegenerative disorders (Alzheimer’s and Parkinson’s disease) and several types of cancer (Bhoola et al. 1992; Jaffa et al. 1992; Ogawa et al. 2000; Shimizu-Okabe et al. 2001; Sharma 2003; Diamandis et al. 2004; Kontos and Scorilas 2012; Fischer and Meyer-Hoffert 2013; Fuhrman-Luck et al. 2014; Prassas et al. 2015). In this scope, KLK3, also known as the prostate-specific antigen (PSA), is by far the most studied member of the KLK family, given its relevance as a tumor marker for prostate cancer. However, recent studies have shown that other KLKs have been proposed as potential biomarkers mainly because of their deregulation in different types of cancer (Clements et al. 2004; Kontos and Scorilas 2012; Fuhrman-Luck et al. 2014). Chapter 1 | General Introduction 19 Figure 4 –KLKs expression patterns in adult tissues. mRNA concentration for each KLK (as indicated in top row) and tissue. The color code at the bottom shows the levels of expression (from Shaw and Diamandis 2007). 2.3.1. Functions in reproductive biology Semen, also known as seminal fluid, is a complex organic medium produced by male reproductive organs and sex accessory glands. It is composed by a wide range of substances that are essential for proper spermatozoa function and fundamental for semen coagulation and liquefaction processes (Jequier 2000). The epididymal fluid that contains the spermatozoa represents only a small fraction of the total ejaculated volume (<5%), whereas the seminal vesicles secretions, including major structural proteins, account for the largest semen fraction (approximately 65%). The remaining portion of semen volume comes mostly from prostate secretions, in a fluid enriched in proteases and Zn2+ Chapter 1 | General Introduction 20 (approximately 30%) (Lundwall and Brattsand 2008). At the moment of ejaculation, in normal physiological conditions, prostate and seminal vesicle secretions are mixed with the epididymal fluid to form a coagulum in the vaginal cavity (Michael et al. 2006; Malm et al. 2007). This gelatinous mass results from the cross-linking of semenogelin 1 (SEMG1), semenogelin 2 (SEMG2) and fibronectin (FN), which are the predominant structural proteins of the seminal plasma, and has the role of entrapping and protecting the spermatozoa. Later, the liquefaction of the coagulum (5-20 minutes after ejaculation) allows a progressive release of the motile spermatozoa together with small peptides that protect sperm cells with their antibacterial, antiviral and antifungal properties (Lilja and Laurell 1985; Lilja et al. 1989; Malm et al. 1996; Peter et al. 1998; de Lamirande et al. 2001; Michael et al. 2006; Edstrom et al. 2008; Lundwall and Brattsand 2008; Zhao et al. 2008). The liquefaction process (Figure 5), involving a stepwise cleavage of SEMG1, SEMG2 and FN, is essentially driven by KLK3 and KLK2 (Lilja 1985; Lilja et al. 1987; Deperthes et al. 1996). In the prostate, the high concentration of Zn2+ regulates the KLK activity with a reversible and allosteric inhibition (Jonsson et al. 2005). Upon ejaculation, the Zn2+ is redistributed through the SEMGs, which consequently activates KLK3 and KLK2 allowing the liquefaction of the semen coagulum. However, as SEMGs are hydrolyzed, Zn2+ is gradually released and the KLKs activity is downregulated again. This mechanism of negative feedback, controlled by the availability of Zn2+, is fundamental to prevent excessive proteolysis which may damage spermatozoa integrity (Robert and Gagnon 1999; Emami and Diamandis 2007). At this stage, a number of endogenous inhibitors and regulatory feedback loops tightly control the cleaving process. For instance, protein C inhibitor (SERPINA5) and eppin (SPINLW1) are two serine protease inhibitors that are known to form complexes with SEMGs, preventing a premature proteolytic cleavage of the semen coagulum by KLK3 (Suzuki et al. 2007; Wang et al. 2007; McCrudden et al. 2008). Importantly, recent studies are proposing that other members of the KLK family beside KLK3 might also play a role in the cascade of semen liquefaction. For example, KLK2, KLK4, KLK5, KLK14 and KLK15 were shown to activate pro-KLK3 in vitro (Deperthes et al. 1996; Takayama et al. 2001a; Takayama et al. 2001b; Michael et al. 2006; Emami and Diamandis 2008) and KLK5 and KLK14 were also reported to cleave SEMGs and FN in vitro (Michael et al. 2005; Michael et al. 2006; Emami et al. 2008). Chapter 1 | General Introduction 21 Despite the recent progress in the understanding of KLKs physiological functions in semen liquefaction, their pathological relevance in male infertility remains largely unknown. Currently, it is estimated that in Western countries approximately 15-20% of couples within reproductive age experience difficulties in achieving pregnancy after one year of regular sexual intercourse. Furthermore, among these couples the male factor Figure 5 - Schematic representation of the semen liquefaction proteolytic cascade. (A) In normal physiologic conditions, KLKs are activated in the prostate through a zymogen activation cascade. KLK activation by other KLK is represented by straight arrows and auto-activation ability is illustrated by curved arrows. The pro-peptide is represented by the yellow rectangle. (B) Upon ejaculation, the sperm-rich epididymal fluid is mixed with prostatic fluids (including KLKs) along with secretions of the seminal vesicles (including SEMG1, SEMG2 and FN), forming the semen coagulum. SEMGs chelate Zn2+ ions, which leads to KLK reactivation and subsequent proteolysis of the SEMGs and FN, resulting in seminal coagulum liquefaction (adapted from Michael et al. 2006 and Prassas et al. 2015). Chapter 1 | General Introduction 28 KLK2 mutation affecting the catalytic triad (D102A) and, consequently, enzyme activity was proposed to contribute to the different semen physiology of this species, in which the semen does not liquefy but forms a copulatory plug instead (Clark and Swanson 2005). Given that the rhesus monkey is a polygamous species in nature, the presence of a copulatory plug is important for sperm competition and mate guarding. However, taking into account that SEMG1 has been inactivated by a frameshift mutation in this species too, KLK2 loss-of-function could result in a reduced ability to dissolve the semen coagulum, thus allowing the formation of a copulatory plug. Still, the role of KLKs in the reproductive system and the dynamical implications of mating behavior remain poorly characterized. Chapter 2 Aims Chapter 2 – Aims 31 The importance of KLKs in the cascade of semen liquefaction together with the evidence that proteolytic and reproductive genes might have been targeted by natural selection, at interand intraspecific level, has motivated the characterization of KLK sequence variation in different primate species, and in healthy and diseased human populations. Specific aims of this work: 1. Unravel the evolutionary history of the most recent KLK duplicates, KLK2 and KLK3, and address a possible correlation to primate mating systems and sperm competition. In a systematic analysis of KLK2 and KLK3, using comparative and phylogenetic procedures, a total of 22 primate species with diverse mating systems and different patterns of semen coagulation were investigated. 2. Characterize a signature of natural selection shaping KLK cluster diversity in Asian populations. The hypothesis of a potential signature of natural selection among KLK genes in Asians, as previously highlighted by independent genome-wide scans of positive selection, was accomplished by undertaking a comprehensive survey of 1000 Genomes (phase I) data and combining it with in vitro functional assays for the most likely candidate variants. Furthermore, a possible ascertainment bias of the 1000 Genomes data was evaluated by Sanger sequencing of several genomic segments across the KLK cluster in a subsample of individuals screened by the 1000 Genomes project. 3. Assess the impact of KLK sequence variation in different infertility phenotypes. An association study of male infertility centered in the genetic screening of the KLK cluster, along with their targets, SEMG1 and SEMG2, and their potential inhibitors of the whey acidic protein four-disulfide core domain (WFDC) locus was performed in a cohort of Portuguese infertility cases and controls. The survey for potential candidate variants for male infertility was achieved through a combination of pooled sample high-throughput sequencing, Sanger sequencing and other genotype screening methods. Chapter 3 Papers Paper I - Birth-and-Death of KLK3 and KLK2 in Primates: Evolution Driven by Reproductive Biology Genome Biol. Evol. 2012 4(12):1331-1338 Chapter 3 | Papers 37 Chapter 3 | Papers 44 Paper II - Adaptive Evolution Favoring KLK4 Downregulation in East Asians Mol. Biol. Evol. 2015 Epub ahead of print Chapter 3 | Papers 47 Chapter 3 | Papers 48 Chapter 3 | Papers 49 Chapter 3 | Papers 50 Chapter 3 | Papers 51 Chapter 3 | Papers 52 Chapter 3 | Papers 53 Chapter 3 | Papers 60 Chapter 3 | Papers 61 Chapter 3 | Papers 62 Paper III - Rare and common variants in KLK and WFDC gene families and their implications into semen hyperviscosity and other male infertility phenotypes In preparation Chapter 3 | Papers 65 Rare and common variants in KLK and WFDC gene families and their implications into semen hyperviscosity and other male infertility phenotypes Patrícia Isabel Marques1,2,3,4, Filipa Fonseca1,2, Ana Sofia Carvalho5, Diana A. Puente3, Isabel Damião6, Vasco Almeida6,7, Nuno Barros8, Alberto Barros8,9, Filipa Carvalho9, Rune Matthiesen5, Victor Quesada3, Susana Seixas1,2 1 - Instituto de Investigação e Inovação em Saúde, Universidade do Porto (I3S), Porto, Portugal; 2 - Institute of Molecular Pathology and Immunology of the University of Porto (IPATIMUP), Porto, Portugal; 3 - Department of Biochemistry and Molecular Biology-IUOPA, University of Oviedo, Oviedo, Spain; 4 - Institute of Biomedical Sciences Abel Salazar (ICBAS), University of Porto, Porto, Portugal; 5 – Human Genetics Department, National Institute of Health Dr Ricardo Jorge (INSA), Lisboa, Portugal; 6 – Center of Infertility and Sterility Studies (CEIE), Porto, Portugal; 7 – Department of Biology, Faculty of Sciences, University of Porto, Porto, Portugal; 8 – Center for Reproductive Genetics Alberto Barros, Porto, Portugal; 9 – Department of Genetics, Faculty of Medicine, University of Porto, Porto, Portugal. Corresponding author: Susana Seixas IPATIMUP, Rua Júlio Amaral de Carvalho 45, 4200-135 Porto, Portugal Phone: +351 225570700; Fax: +351 225570799 E-mail: sseix[email protected] Chapter 3 | Papers 66 Abstract The human kallikrein (KLK) and whey acidic protein four-disulfide core domain (WFDC) gene families, located at chromosomes 19q13.3-13.4 and 20q13, respectively, encode molecules with key roles in the cascade of semen coagulation and liquefaction. Semenogelins 1 and 2 (SEMGs), the main components of the semen coagulum, establish a cross-linked matrix that entraps spermatozoa. In contrast, KLK3 and KLK2, further assisted by other KLKs, hydrolyze the coagulum in a process crucial for spermatozoa motility. Here, we examined the contribution of KLK and WFDC variation into human infertility, by performing a pilot screening of coding and non-coding regions, covering approximately 93 kb of genomic sequence, by means of a pooled sample high-throughput approach, later, followed by a genotyping survey for most promising candidates in a cohort of cases (N=238) and controls and (N=217). Among the 456 variants identified in the pooled sequencing, 296 were low-frequency for which a higher burden of deleterious alleles was detected in cases for KLKs. Eleven variants were confirmed to be overrepresented in cases and likely affecting KLK4 and KLK12 structure, KLK3 activity or KLK3, KLK14 and KLK15 gene expression. In SEMGs we identified 3 rare variants expected to modify the profile of cleaved peptides with potential effects in spermatozoa interactions. A common nucleotide substitution in KLK7 (rs1654526) and a copy number variation in SEMG1, were on the other hand associated to a reduced risk for different infertility phenotypes. Overall, the results support the importance of KLKs and SEMGs in male reproduction and provide evidence for a contribution of their genetic variation into semen hyperviscosity and asthenozoospermia. Chapter 3 | Papers 67 Introduction Infertility is a major reproductive disorder characterized as the inability to achieve a viable pregnancy after one year of regular sexual intercourse, which affects up to 20% of couples within reproductive age in Western countries (Rowe et al. 1993; Cedenho 2007). In the disease pathogenesis, male and females are thought to have nearly equivalent contributions, as many physiological processes are required to achieve a successful fertilization. Indeed, in most couples, the two partners are often subfertile and only in a smaller percentage of cases the male factor is considered as the primary cause for a reproductive failure (Thonneau et al. 1991; Practice Committee of American Society for Reproductive 2012). So far, the recognized causes for male infertility are congenital defects leading to reproductive malformations, endocrine and immunologic dysfunction, mechanic trauma, urogenital infections causing post-testicular obstruction, impaired spermatogenesis and also major and minor genetic abnormalities, such as aneuploidies, translocations and deletions in Y and other chromosomes, as well as point mutations in genes like CFTR, AR, DAZL, and SRY (Jungwirth et al. 2012; Tahmasbpour et al. 2014). Although, infertility has been proposed as a complex disease with a strong genetic component, it continues mostly unexplained, as the large proportion of disease cases are frequently defined as idiopathic or have an unknown etiological cause (Carrell and Aston 2011; Aston 2014). In the clinical practice, abnormal semen parameters as determined by reference threshold values of the World Health Organization (WHO), are regarded as evidence for male infertility and accordingly, infertility patients are classified into different non-mutually exclusive phenotypes. These include hyperviscosity, for a persistence of semen viscosity; oligozoospermia or azoospermia, for lower spermatozoa counts; asthenozoospermia, for reduced spermatozoa motility and teratozoospermia, for altered spermatozoa morphology (WHO 1999; WHO 2010). Semen is a body fluid that results from an assorted mixture of the spermatozoa-rich secretions of testis and epididymis, with the products from the seminal vesicles, prostate, and bulbourethral glands. Upon ejaculation, the semen forms a gelatinous mass – the semen coagulum – that entraps the spermatozoa in a matrix of linked seminal proteins. Over a short period of time (5 to 20 minutes), the semen coagulum starts to be liquefied through the activity of several enzymes, allowing a progressive release of the spermatozoa and the regain of their motility (Lilja 1985; Deperthes et al. 1996; Robert et al. 1997). Chapter 3 | Papers 68 Several proteins known to play key roles in the cascade of semen coagulation and liquefaction belong to the kallikrein (KLK) gene family, found in chromosome 19q13.3-13.4 region, and to the whey acidic protein four-disulfide core domain (WFDC) family, located at chromosome 20q13. Specifically, the two major structural seminal proteins enrolled in the semen coagulation and spermatozoa immobilization, semenogelin 1 (SEMG1) and semenogelin 2 (SEMG2), belong both to the WFDC locus (Lundwall and Brattsand 2008; Clauss et al. 2011). On the other hand, the core enzymes involved in the semen liquefaction and SEMGs hydrolysis into multiple shorter peptides are KLK3 (also known as prostate-specific antigen – PSA) and KLK2 (Lilja 1985; Lilja et al. 1987; Deperthes et al. 1996). Moreover, the WFDC cluster also encompasses 17 small serine protease inhibitor genes, which include SLPI and PI3, two ubiquitous molecules with antimicrobial activities at reproductive mucosal surfaces, and EPPIN, a molecule that coats the spermatozoa in a protein complex with fibronectin and SEMGs, modulates KLK3 activity and coagulum proteolysis, and also protects spermatozoa form bacterial attacks (Yenugu et al. 2004; Wang et al. 2005; Williams et al. 2006; Wang et al. 2007; Weldon et al. 2007; McCrudden et al. 2008; Zhao et al. 2008; Zhang et al. 2013). Less is known about other WFDC molecules, but most genes were found to be mainly expressed in male reproductive tissues and several were detected in the human seminal plasma by proteomic profiling (Clauss et al. 2002; Thimon et al. 2008; Batruch et al. 2012; Chhikara et al. 2012). Conversely, the KLK locus includes a total of 15 trypsinor chymotrypsin-like serine protease genes (KLK1-KLK15) with pervasive activities in diverse proteolytic cascades, including semen liquefaction. As example, KLK4, KLK5, KLK14 and KLK15 were all shown to regulate the activity of KLK3, and KLK5 and KLK14 were reported to overlap with KLK3 in the hydrolysis SEMGs and fibronectin (Takayama et al. 2001a; Takayama et al. 2001b; Michael et al. 2006; Emami et al. 2008; Emami and Diamandis 2008). In the latest years, the evidence supporting a possible contribution of KLK and WFDC families into different male infertility phenotypes has been consolidating. In the study by Emami et al., most KLKs were correlated to a downregulation of protein expression and a hyperviscosity phenotype, whereas the delay of semen liquefaction was associated to a restricted set of proteins (KLK2-3 and KLK13-14) displaying reduced protein expression (Emami et al. 2009). Furthermore, in the same study a link between KLK14 expression levels and asthezoospermia was established, and more recently KLK3 was found to be significantly upregulated in patients combining oligo and teratozoospermia phenotypes (Emami et al. 2009; Sharma et al. 2013). On the other hand, SEMG1 has been recently shown to be upregulated in infertile patients with and without asthenozoospermia, but no correlation could be observed between the patterns of Chapter 3 | Papers 69 SEMG degradation and hyperviscosity (Martinez-Heredia et al. 2008; Esfandiari et al. 2014; Legare et al. 2014; Yu et al. 2014). Moreover, in the EPPIN gene two single nucleotide variants (SNV) were shown to correlate with semen quality in the Han-Chinese population, one presenting a decreased risk to low sperm number and the other an increased risk to abnormal motility (Ding et al. 2010a; Ding et al. 2010b). In this study, we sought to investigate in which extent the genetic variation within KLK and WFDC families affects the regular process of semen coagulum liquefaction, and underlies hyperviscosity and other infertility phenotypes. By performing a comprehensive survey of KLK and WFDC coding and non-coding regions using a high-throughput sequencing strategy, we demonstrated that KLKs have an excess of low-frequency variants among infertility cases and we validated a total of 12 candidate variants of male infertility in KLKs but also in WFDCs, with expected impact in the proteolytic processing of the semen coagulum. Chapter 3 | Papers 76 inference) nature of the method implemented for the screening of KLK and WFDC gene variability. Globally, in this phase of the study (phase I), a total of 456 SNVs were identified (Supplementary Table S2), in which 296 (64.9%) were low-frequency (MAF < 0.05) and 104 (22.8%) were novel variants. In addition, 98 (21.5%) SNVs were located in coding exons (58 nonsynonymous and 40 synonymous), 72 (15.8%) in UTRs and 16 (3.5%) in splice regions. Analysis of common variants association to male infertility To test the association to male infertility of the 160 identified common SNVs (MAF  0.05), we performed a series of comparisons between controls and different groups of cases: all infertility phenotypes (HV+NV), hyperviscosity cases (HV), or asthenoand oligozoopermic phenotypes without hyperviscosity (NV). Precisely, 29 SNVs, 23 in KLKs and 6 in WFDCs, were found to have significant associations in at least one of the comparisons (Supplementary Table S3). However, only 2 variants in the HV group maintained their statistical association after controlling for multiple test comparisons (P < 0.0003125). The SNV showing the strongest value (P = 0.0002) was a synonymous substitution in WFDC6 (rs41304411, ENST00000372665.3:c.366C>A, p.I122I), with no obvious functional effect. The other associated SNV, rs1654526 (ENST00000595820.5:c.606+585C>T) (P = 0.0003), was located in KLK7, within a enhancer region that harbors several binding sites for transcription factors (FOS, FOSL2 and JUND) and repressors (CTCF) as indicated by chromatin segmentation and Chip-seq data from ENCODE, respectively (Supplementary Fig. S5). Interestingly, both SNVs showed higher frequencies in controls than in cases, which suggest an undirected link or a possible protective role of these SNVs to male infertility. Burden tests of low-frequency variants The vast majority of SNVs identified in our study were low-frequency variants (MAF < 0.05), for which standard statistic tests have low power to detect associations, especially in relative small sample sizes (few hundred individuals). To circumvent this limitation and to investigate whether there is a higher burden of low-frequency deleterious variants (nonsynonymous and splice region SNVs) in cases than in controls, we applied the Calpha statistic under several gene combined analysis (Neale et al. 2011) (Table 1). We started by analyzing KLK and WFDC genes altogether in a single group, for which a significant enrichment of low-frequency variants in infertility cases was detected independently of the considered disease phenotype (HV, NV and HV+NV). In two other Chapter 3 | Papers 77 tests, to address if KLK and WFDC clusters were differentially enriched in low-frequency SNVs, we calculated the statistic for each gene family, separately. This analysis confirmed a higher burden of low-frequency variants for all phenotypes, but solely for KLK genes. An equivalent approach was used to inquiry if cases had a higher burden of regulatory SNVs (5’ and 3’ UTR variants) than controls, nevertheless, none of the analysis yielded significant results (Table 1). Candidate variant genotyping To confirm variant calling and allele frequency estimates obtained in the first phase of the study (phase I), we performed a genotyping screening by Sanger sequencing for the most promising candidate genes of male infertility (phase II; Supplementary Fig. S1). Here, we prioritized the analysis of several gene regions based in the following criteria: 1) previous evidence of deregulated activity in infertility cases and/or recognized role in male reproduction; 2) gene contains at least 2 low-frequency variants with potential deleterious effects; and 3) SNV is absent in controls or displays major frequency differences between cases and controls (two or more times higher). Specifically, we selected for the Sanger sequencing study segments of KLK3, KLK4, KLK6, KLK12 and KLK14 genes. The function of KLK12 has not been fully elucidated yet, but according to its mRNA expression in male genital tissues and protein identification in seminal plasma, it is likely to have an important function in male reproductive biology (Shaw and Diamandis 2007). As a whole, the results of the Sanger sequencing screening provided a good fit to the estimated MAFs in pooled samples, as demonstrated by the strong correlations between SNVs for the three sample groups considered (HV: r2 = 0.9725; NV: r2 = 0.9695; controls: r2 = 0.9509) (Fig. 2). Noticeably, in these segments the estimates of allele frequency derived from the pooled sequencing had similar levels of accuracy for both common (MAF ≥ 0.05) and low-frequency (MAF < 0.05) variants (Fig. 2). In the category of low-frequency variants, several SNVs emerged as promising candidates for male infertility (Table 2 and Supplementary Fig. S6). In KLK3, we confirmed for the HV group the segregation of rs111901464 (ENST00000593997.5: c.658G>A, p.E220K), a SNV located in the intron IV of KLK3 also predicted to have a negative impact in a shorter KLK3 isoform; an increased incidence in HV cases of a nonsynonymous substitution (rs61729813, ENST00000326003.6: c.629C>T, p.S210W); and we identified in a single asthenozoospermic patient, a variant previously uncovered by the pooled sequencing (rs182759459, ENST00000326003.6: c.391G>A, p.E131K). While both p.E220K and p.E131K entailed a substitution of two amphipathic residues with opposite side chain charges that are frequently involved in salt-bridges and in interactions Chapter 3 | Papers 78 with non-protein atoms, the p.S210W involved the change of a small polar amino acid to a large hydrophobic aromatic residue. In KLK14, we confirmed an increased prevalence among infertility cases of 2 deleterious SNVs, a splice donor site (rs117229324, ENST00000391802.1: c.26+1G>A) likely affecting normal splice processing (HSF and MaxEnt scores changes from 93.16 to 66.33 and 10.24 to 2.06, respectively) and a nonsynonymous variant (rs112658494, ENST00000391802.1: c.412C>T, p.R138W) replacing a positively charged amino acid by a large hydrophobic aromatic residue rarely engaged in non-protein atom binding. Furthermore, in KLK12, we found 2 SNVs confined to infertility cases, one identified in a single patient and causing a loss of a disulfide bond (rs140609488, ENST00000319590.8:c.587G>A, p.C196Y) and the other detected in a few cases leading to a substitution of an extremely conserved proline (rs61742847, ENST00000319590.8: c.101C>T, p.P34L) across KLK family (Supplementary Fig. S6). Finally, in KLK4 and KLK6, we validated the presence of two rare deleterious variants found in isolated cases. Whereas, for KLK4, we uncovered in an asthenozoospermia case a novel p.Q42L replacement affecting a highly conserved amino acid among KLKs (Supplementary Fig. S6), in KLK6, we disclosed for a combined phenotype of asthenozoospermia and hyperviscosity a p.T234M (rs77760094, ENST00000310157.6: c.701C>T) mutation located in the C-terminal region. Noticeably, for KLK6 we also found in two controls a novel variant p.I216N with predicted functional consequences. The p.I216N replaces a hydrophobic amino acid with an aliphatic side chain by a polar residue frequently involved in protein activity and binding sites, such as the catalytic triad of several cysteine proteases. To investigate the putative structural consequences of the candidate variants, we mapped the affected residues onto the solved three-dimensional (3D) structure of KLK4 bound to a competitive inhibitor. Notably, the alignment of kallikrein sequences showed that several variants cluster at close or even equivalent positions in the primary structure (Fig. 3A). Moreover, the threading of the alignment onto a 3D structure revealed further clustering of positions 196 of KLK12 and 210 of KLK3 to position 215 in KLK6 (Fig. 3B). Interestingly, this cluster is located in the inhibitor-binding pocket of the structure, which suggests that these variants may directly affect the binding of the corresponding kallikreins to their substrates. In fact, the C196Y variant in KLK12 is expected to destroy a disulfide bond, termed SS6, necessary for the catalytic activity of kallikreins (Oka et al. 2002). The remaining variants map to different parts of the structure, and still show a tendency to cluster. Thus, variants E131K in KLK3 and R138W in KLK14 map to the same position in the structure (Fig. 3B). In the absence of further functional information, this suggests that variants in this position are selected for in infertile patients. Similarly, Chapter 3 | Papers 79 variants Q42L in KLK4 and P34L in KLK12 cluster in consecutive positions in the structure. Both positions are well conserved among kallikreins (Fig. 3A and Supplementary Fig. S6) and appear about ten residues downstream from the activation site, suggesting that these variants might affect protein folding and the activation process. Finally, the T234M variant in KLK6 is located in a -helix close to the C-terminus of the protein. Beside the validation of candidate genes by Sanger sequencing 3 more dispersed variants were selected for genotyping using a SNaPShot multiplex reaction approach. These included a SNV located in the KLK8 5’UTR (rs74705037, ENST00000600767.5: c.- 29C>T), within a region recognized by ENCODE chromatin segmentation data as an enhancer and predicted by MatInspector to contain several binding motifs for Zinc finger proteins (loss of Zinc finger and BTB domain-containing protein 7A and Zinc finger protein GLIS2 binding sites). Another SNV allocated to a KLK15 splice region (rs3212852, ENST00000598239.5: c.481+5G>A) probably impairing normal mRNA processing (HSF and MaxEnt scores changes from 89.08 to 76.92 and 9.49 to 4.62, respectively). A latter SNV placed in EPPIN intron I (rs75681320, ENST00000354280.8: c.92-438C>T) in an insulator region containing several transcription binding factors (CTCF, SMC3, RAD21) as identified by ENCODE chromatin segmentation and Chip-seq data, respectively. Extended association study of male infertility To increase the statistic power of our study, we carried out an extended genotype analysis for the most promising candidate variants for male infertility. In this phase (phase III), we screened an additional panel of 95 cases (36 HV and 59 NV) and 138 controls (34 individuals with normal semen parameters, 10 fertile men and 94 random Portuguese males), for 7 low-frequency SNVs and a common SNV (Supplementary Fig. S1). The selected low-frequency SNVs included: 3 nonsynonymous substitutions confirmed by Sanger sequencing and predicted as deleterious (KLK3 p.E220K - rs111901464; KLK12 p.P34L - rs61742847 and KLK14 p.R138W - rs112658494; the splice donor site (KLK14: c.26+1G>A - rs117229324); and the 3 SNVs evaluated by SNaPshot (KLK8: rs74705037; EPPIN: rs75681320 and KLK15: rs3212852). The single common variant selected for the extended screening was the one located on KLK7 (rs1654526). The case-control analysis of the entire dataset comprising a total of 238 infertility cases (111 HV and 127 NV) and 217 controls (Table 3 and Supplementary Fig. S1) corroborated a trend toward increased MAFs in infertility cases for all low-frequency variants (Table 3). However, the sample size was not enough to reach statistical significance in most circumstances. Indeed, only the KLK12 p.P34L (rs61742847) Chapter 3 | Papers 80 replacement showed a significant association for both HV+NV and HV group comparisons (P = 0.0388 and P = 0.0384, respectively), thus suggesting a possible contribution into the hyperviscosity phenotype. Nonetheless, the candidate variants identified in KLK3, KLK15 and EPPIN, were found to present at least three times higher frequencies in the HV group than in controls. In contrary, the KLK14 splice variant showed an equivalent frequency increment in the NV group that could be connected to an asthenozoospermia phenotype, including in the single case identified among the HV group. Finally, the common KLK7 (rs1654526) was confirmed as significantly associated with a reduced susceptibility to semen hyperviscosity (P = 0.0035) and male infertility (P =0.0258). SEMGs genetic screening The survey of SEMGs was centered in exon II, which covers nearly all protein coding sequence, except for a short N-region included in exon I. Similarly to the approach used in the high-throughput sequencing study, we started by analyzing the same cohort of 143 cases and 79 controls by Sanger sequencing. In the first analysis of SEMGs variation, we identified 6 and 4 SNVs in SEMG1 and SEMG2, respectively (Supplementary Table S4). Among the SNVs found in SEMGs, only rs147894843 (ENST00000372781.3: c.1199G>A, p.G400D) in SEMG1, and rs2233903 (ENST00000372769.3: c.835C>T, p.H279Y), rs2071650 (ENST00000372769.3: c.1102G>C, p.G368R) and rs139977707 (ENST00000372769.3: c.1654G>C, p.E552Q) in SEMG2 were predicted to affect protein function. In addition, in the SEMG1 we also covered the previously described CNV corresponding to the 5 or 6 repeat units (Jensen-Seaman and Li 2003; Lundwall et al. 2003; Miyano et al. 2003). To investigate whether there was a higher burden of deleterious SNVs among SEMGs in cases than in controls, we applied again the C-alpha statistic (Neale et al. 2011) (Supplementary Table S5). However, we did not observe any variant enrichment in all tests performed. In this initial phase, a single silent substitution located in SEMG1 and restricted to controls (p.T293T rs17850164) was found to displayed a significant association (P = 0.0423, in the HV+NV comparison). Nevertheless, 2 nonsynonymous substitutions were absent in controls (SEMG1 p.G400D and SEMG2 p.E552Q); 2 linked variants in SEMG2 (p.H279Yand p.G368R) were slightly increased in cases, and the 5 repeat allele of SEMG1 showed a higher frequency in controls than in cases. Interestingly, the significantly associated p.T293T mutation had no predicted effect on mRNA secondary structure but it appeared as linked to the 5-repeat allele. Thus, in the extended study, we decided to genotype for SEMG1 the p.G400D replacement and the CNV; and for SEMG2 the p.H279Y and p.E552Q substitutions. Chapter 3 | Papers 81 The analysis of the full case-controls datasets revealed a significant association of SEMG1 p.G400D (rs147894843) variant in HV+NV and NV groups (P = 0.0388 and P = 0.04920, respectively) (Table 4). This variant creates a potential cleavage site (P1-P1’: DE search in MEROPS) for cathepsin D and metalloprotease 2, two proteases found at high and low abundances in semen, respectively and also for caspase-3, a cysteine protease previously associated to male infertility and asthenozoospermia (Almeida et al. 2005; Rawlings et al. 2014). Consistently, in our study 4 out of 5 cases carrying the p.G400D had an asthenozoospermia phenotype, for which a significant association was also obtained (P = 0.0442, asthenozoospermia vs. controls). The SEMG2 variants did not reach significance in our group comparisons, nevertheless, we found a slighter increment on the p.H279Y (rs2233903) replacement in the NV group and a 7 times augmented frequency of p.E552Q (rs139977707) in the HV. Oddly, this SNV is likely to generate a novel cleavage site for KLK3 (P1-P1’: QS) in the SEMG2 sequence (Malm et al. 2000). The CNV maintained the tendency toward a lower frequency of the 5-repeat allele in infertile patients (P = 0.0667 for HV+NV cases vs. controls), as previously reported in Asians (Miyano et al. 2003). Notably, in our sample, we could detect a significant lowering of the 5-repeat allele in oligozoospermia patients (P = 0.02928) but not in asthenozoospermia as it has been previously hypothesized (Miyano et al. 2003). At last, in the extended survey of SEMGs we also discovered in a single individual exhibiting a combined phenotype of hyperviscosity and asthenozoospermia a novel SEMG1 variant (hg19 chr20: g.43837061T>C, p.Y315H), expected to abolish one of the KLK3 cleavage sites (Rawlings et al. 2014). Proteomic validation To assess a possible deleterious effect of the candidate variants that could result in the degradation of the mutant protein, we carried out a proteomic screening of KLK3 p.E131K and p.S210W substitutions in the seminal plasma of individuals bearing these mutations. In both cases, it was possible to identify the mutant variant in comparable levels to those for the wild type (Fig. 4 for p.S210W; for p.E131K data is not shown), indicating that the mutant allele is secreted and not degraded by the cellular machinery. However, other deleterious functional effect of these variants in KLK3 function cannot be excluded: the p.S210W substitution by being placed in the binding pocket can still affect the interactions to substrates, and in the p.E131K the loss of glutamate residue could potentially affect the binding to ions like zinc, an important modulator of KLK activity in the seminal plasma. Chapter 3 | Papers 82 Discussion We performed a comprehensive study of male infertility focused on the genetic variation of KLK and WFDC clusters, by applying a next-generation sequencing approach to the analysis of pooled DNA samples, allowing in a cost-effective manner to identify multiple susceptibility markers in cases with and without hyperviscosity. A follow-up genotype screening of selected locus regions confirmed the accuracy of the method, which permitted the variant calling and frequency inference of both common and rare alleles. Globally, an enrichment of potential candidate variants (high and low-frequency) was detected in KLKs in comparison to WFDCs, supporting a greater impact of the first cluster in human reproduction and fertility. Furthermore, in agreement with KLKs higher levels of sequence variation in infertility cases, other authors had reported a consistent downregulation of protein expression for most tested KLKs (KLK1-3, KLK5-10 and KLK1314) in the semen of individuals with abnormal viscosity and liquefaction parameters (Emami et al. 2009). Although, a greater emphasis was given in our study, in particular in the genotyping surveys, to nonsynonymous and splice variants, these may only explain a small fraction of KLK quantitative differences in semen through mechanisms of abnormal protein and mRNA degradation. Indeed, among validated variants only the p.P34L substitution in KLK12, which is located in a highly conserved residue across the entire family, as well as, in other far related serine proteases, and the mutated donor splice region of KLK15, could hypothetically be linked to hyperviscosity by those straightforward effects. However, neither KLK12 nor KLK15 were previously evaluated in the semen, and whereas the later protease is known to activate pro-KLK3 in-vitro, a role of KLK12 in male reproduction was not yet been elucidated (Takayama et al. 2001a). Worth to note, in this study two common variants in KLK12 region that did not surpass multiple testing were also associated to the hyperviscosity phenotype. On the other hand, SNVs positioned in regulatory regions may also cause significant differences in protein expression, if these disrupt a binding motif for an important transcription factor. Even though, this category of genetic variants of male infertility has been less explored in the genotyping phases of the study because of the less confident nature of bioinformatic predictions, several candidates are likely to have their link to hyperviscosity explained by such phenomena. This is likely to be the case of the low-frequency variant in the 5’ UTR of KLK8 (rs74705037), slightly augmented in hyperviscosity cases and the significantly associated variant, located in an intron of KLK7 Chapter 3 | Papers 83 (rs1654526), both placed in enhancer regions previously shown by ENCODE to interact with multiple regulatory elements. Interestingly, KLK7 besides being found in extremely reduced concentrations in hyperviscosity cases (Emami et al. 2009), its gene harbors three other common SNVs less strongly associated to this phenotype (rs1991820, rs1991819, and rs1991818). Another genotyped SNV possibly correlated to the hyperviscosity phenotype through gene downregulation is a low-frequency variant in KLK3 (rs111901464) positioned in a region identified by ENCODE chromatin segmentation as weaker enhancer, but that could also affect a shorter protein isoform lacking the catalytic triad (p.E220K). Up till now, a reduced expression of KLK3 has only been described in cases of delayed viscosity and displaying reduced spermatozoa number with abnormal morphology (Emami et al. 2009; Sharma et al. 2013). Two other nonsynonymous variants, predicted as deleterious, were found to be slightly increased in hyperviscosity (near two times higher), the KLK3 p.S210W (rs61729813) and the KLK14 p.R138W (rs112658494). However, these SNVs if truly linked to infertility are only expected to cause qualitative changes in protein interactions with other molecules. Specifically, the p.S210W variant of KLK3 was confirmed to not differentially affect the protein content in the semen and since it is located in a balcony region of the catalytic pocket (Ser210-Gly225) exposed to bulk solvents probably affects protease activity (Debela et al. 2006). Here, a substitution of a small polar residue (serine) by a large aromatic amino acid (tryptophan), where other KLK often present a glutamine (Q), is likely to restrain the substrate acceptance in the KLK3 catalytic pocket. Despite the limited number of SNVs identified in the WFDC cluster, two possible susceptibility variants were discovered to be more prevalent among the hyperviscosity cases, a SNV located in a regulatory region, specifically in an insulator of EPPIN, a smaller protease inhibitor reported to target KLK3, also implicated in antimicrobial activities (McCrudden et al. 2008) and a p.E552Q substitution in SEMG2. This later variant is likely to introduce a novel cleave site for KLK3 in the C-terminal region of the protein, but so far no specialized function has been attributed yet to such region (Robert et al. 1997; Robert and Gagnon 1999). Still, other susceptibility variants evaluated in the genotyping screening were overrepresented among infertility cases with normal viscosity. These included a nonsynonymous variant located in a critical structural region, p.Q42L in KLK4, and the disrupted donor splice of KLK14, all possible contributing to a lowering of the KLK content. Conversely, the remaining identified variants are mostly likely affecting molecular Chapter 3 | Papers 84 interactions or protein activity, namely the KLK3 p.E131K, the KLK12 p.C196Y and the SEMG1 substitutions p.Y315H and p.G400D. Notably, both SEMG1 variants are expected to alter the protein proteolytic processing, but only the later substitution is in close proximity to a key sequence recognized as a thyrotropin-releasing hormone (TRH) like peptide (375-397 residues). Despite some controversy about the SEMG1 origin of the TRH-like peptides found in the human semen, these were demonstrated to increase the capacitation of spermatozoa (Khan and Smyth 1993; Huber et al. 1998; Robert and Gagnon 1999). More recently, nearly the same sequence of SEMG1 (376-388 residues) was found to bind to the CD52 glycosylphosphatidylinositol anchored antigen presented by spermatozoa, which also takes part in semen coagulation and its released during liquefaction (Flori et al. 2008). Most of these susceptibility SNVs found in normal viscosity cases were otherwise correlated in several instances to asthenozoospermia, a finding that may be still consistent with abnormal patterns of semen liquefaction. SEMGs and specially SEMG1 are described to bound spermatozoa and to modulate their motility in a doseand time dependent fashion (Robert and Gagnon 1996; Yoshida et al. 2008; Mitra et al. 2010). So far, the motility inhibitory peptides were correlated with N-terminal peptides released during semen liquefaction (-inhibin-92 and -inhibin-31) and containing cysteine 239 residue, but the functional relevance of other SEMG1 regions should not be discarded (Silva et al. 2013). There is a growing body of evidence for the interaction of SEMG1 with other biomolecules, in most cases, with functional outcomes in spermatozoa motility. SEMG1 has been reported as a target for S-nitroso-glutathione, and for prolactin inducible protein, to bind zinc and while in EPPIN-complexes to interact with spermatozoa calcium channels (Lefievre et al. 2007; Yoshida et al. 2008; O'Rand and Widgren 2012; Tomar et al. 2013). Conversely, in asthenozoospermia patients SEMG1 was shown to remain bound to spermatozoa, to be increased in their semen, and its mRNA to be highly expressed by spermatozoa (Zhao et al. 2007; Martinez-Heredia et al. 2008; Terai et al. 2010; Yu et al. 2014). In overview, all identified variants are expected to have a negative impact in fertility through modifications of semen liquefaction process, resulting from a deregulation of protein/peptide activity and concentration levels. Still, the consequences in the proteolytic cascade of each SNV might be difficult to disentangle, given that KLKs are known to have pervasive links with each other, often activating other members and overlapping in substrate affinity (Lawrence et al. 2010). Furthermore, low-frequency variants were always found in single heterozygosity cases, which indicate that in those individuals they are only Chapter 3 | Papers 85 playing a part into the genetic background of male infertility. Moreover, candidate SNVs are only present in small fraction of cases and therefore, far from fully explaining KLK and SEMG association to hyperviscosity and asthenozoospermia phenotypes. This finding is concordant with a view of male infertility as a complex disease resulting from a complex network of genetic and non-genetic factors with a wide range of susceptibility effects. Chapter 3 | Papers 92 Zhao C, Huo R, Wang FQ, Lin M, Zhou ZM, Sha JH. 2007. Identification of several proteins involved in regulation of sperm motility by proteomic analysis. Fertil Steril 87:436-438. Zhao H, Lee WH, Shen JH, Li H, Zhang Y. 2008. Identification of novel semenogelin Iderived antimicrobial peptide from liquefied human seminal plasma. Peptides 29:505-511. Chapter 3 | Papers 93 Table 1 – Burden tests for KLK and WFDC low-frequency variants. Cases (HV+NV) vs. Controls HV Cases vs Controls NV Cases vs. Controls Nonsynonymous and splice region SNVs KLKs and WFDCs C-alpha = 50.9067 (P-value = 0.0327) C-alpha = 30.7813 (P-value = 0.0034) C-alpha = 22.9367 (P-value = 0.0080) WFDCs C-alpha = 9.1536 (P-value = 0.2872) C-alpha = -1.1391 (P-value = 0.8730) C-alpha = 6.5306 (P-value = 0.1536) KLKs C-alpha = 41.7530 (P-value = 0.0106) C-alpha = 31.9205 (P-value = 0.0001) C-alpha = 16.4061 (P-value = 0.0066) UTR SNVs KLKs and WFDCs C-alpha = -11.0989 (P-value = 0.5599) C-alpha = -2.0320 (P-value = 0.8988) C-alpha = -9.2201 (P-value = 0.5252) WFDCs C-alpha = 7.2740 (P-value = 0.0696) C-alpha = 5.0348 (P-value = 0.0654) C-alpha = 3.4013 (P-value = 0.1822) KLKs C-alpha = -18.3729 (P-value = 0.1808) C-alpha = -7.0668 (P-value = 0.6191) C-alpha = -12.6214 (P-value = 0.2766) Significant P-values (P < 0.05) are underlined Chapter 3 | Papers 94 Table 2 – Low-frequency variants surveyed in phase II. Gene SNP ID Genomic position (hg19) Consequence MAF (number of chromosome) Control Cases (HV+NV) HV cases NV cases KLK3 rs182759459 51361469 p.E131K 0.000 (0) 0.003 (1) 0.000 (0) 0.007 (1) KLK3 rs61729813 51361850 p.S210W 0.013 (2) 0.014 (4) 0.020 (3) 0.007 (1) KLK3 rs111901464 51361879 intronic* 0.000 (0) 0.010 (3) 0.020 (3) 0.000 (0) KLK4 rs138071534 51412573 p.E53E 0.000 (0) 0.003 (1) 0.000 (0) 0.007 (1) KLK4 N/A 51412607 p.Q42L 0.000 (0) 0.003 (1) 0.000 (0) 0.007 (1) KLK4 N/A 51412717 intronic 0.000 (0) 0.003 (1) 0.007 (1) 0.000 (0) KLK4 N/A 51412740 intronic 0.000 (0) 0.003 (1) 0.007 (1) 0.000 (0) KLK6 rs77760094 51462454 p.T234M 0.000 (0) 0.003 (1) 0.007 (1) 0.000 (0) KLK6 N/A 51462508 p.I215N 0.013 (2) 0.000 (0) 0.000 (0) 0.000 (0) KLK8 rs74705037 51504808 5’UTR 0.006 (1) 0.035 (10) 0.033 (5) 0.037 (5) KLK12 rs140609488 51534048 p.C196Y 0.000 (0) 0.003 (1) 0.000 (0) 0.007 (1) KLK12 rs61742847 51537332 p.P34L 0.000 (0) 0.010 (3) 0.013 (2) 0.007 (1) KLK14 rs112658494 51582808 p.R138W 0.006 (1) 0.028 (8) 0.040 (6) 0.015 (2) KLK14 rs117229324 51585978 splice donor 0.006 (1) 0.017 (5) 0.007 (1) 0.029 (4) KLK15 rs3212852 51330129 splice region 0.013 (2) 0.017 (5) 0.027 (4) 0.007 (1) EPPIN rs75681320 44174847 intronic 0.013 (2) 0.024 (7) 0.040 (6) 0.007 (1) * p.E220K predicted as damaging in a shorter isoform. N/A – Not applicable, novel variants without ID. Chapter 3 | Papers 95 Table 3 – Combined case-control association analysis from phase III. NV cases P-value 0.6001 0.3019 0.4120 0.3673 0.5653 0.1362 0.2284 0.5077 Significant nominal P-values (P < 0.05) are underlined. MAF 0.004 0.246 0.024 0.004 0.012 0.016 0.016 0.012 HV cases P-value 0.1150 0.0035 0.3236 0.0384 0.1860 0.7343 0.0919 0.0803 MAF 0.014 0.171 0.027 0.014 0.027 0.005 0.023 0.027 Cases (HV+NV) P-value 0.2170 0.0258 0.3207 0.0388 0.3686 0.2653 0.0965 0.1713 MAF 0.008 0.210 0.025 0.011 0.019 0.011 0.019 0.019 Control MAF 0.002 0.267 0.018 0.000 0.014 0.005 0.007 0.009 Consequence Intronic Intronic 5’UTR p.P34L p.R138W splice donor splice region Intronic SNP ID rs111901464 rs1654526 rs74705037 rs61742847 rs112658494 rs117229324 rs3212852 rs75681320 Gene KLK3 KLK7 KLK8 KLK12 KLK14 KLK14 KLK15 EPPIN Chapter 3 | Papers 96 Table 4 – Combined case-control association analysis for SNVs in SEMG1 and SEMG2. NV cases P-value 0.0492 0.1313 0.1581 0.6327 Significant nominal P-values (P < 0.05) are underlined. MAF 0.012 0.028 0.028 0.000 HV cases P-value 0.1142 0.1355 0.6383 0.1150 MAF 0.009 0.027 0.014 0.014 Cases (HV+NV) P-value 0.0388 0.0667 0.2856 0.3475 MAF 0.011 0.032 0.021 0.006 Control MAF 0.000 0.047 0.014 0.002 Consequence p.G400D del Ia (aa320-379) p.H279Y p.E552Q SNP ID rs147894843 5 repeat units rs2233903 rs139977707 Gene SMEG1 SMEG1 SMEG2 SMEG2 Chapter 3 | Papers 97 Figure 1 – Minor allele frequencies (MAFs) from 1000 Genomes data vs. controls from pooled sequencing. Allele frequency estimates for 277 SNVs based on pooled sequencing from the control group were compared with the described European average frequencies from 1000 Genomes project phase III samples. r2 - correlation coefficient (r2 = 0.826). Figure 2 – Minor allele frequencies (MAFs) from pooled sequencing vs. Sanger sequencing. Estimated MAFs based on pooled sequencing is plotted against the actual frequencies as determined by individual Sanger sequencing for the surveyed regions. r2 - correlation coefficients (HV: r2 = 0.9725; NV: r2 = 0.9695; controls: r2 = 0.9509). The data from HV cases, NV cases and controls are represented in orange, blue and green, respectively. Chapter 3 | Papers 98 Figure 3 – Structural characterization of the KLK low-frequency variants. (A) Alignment of the amino acid sequences of the variant kallikreins. Variant sites are framed in red. Complete conservation is shown in dark blue background, whereas partial conservation is shown on a light blue background. The catalytic serine is highlighted with a red arrow. (B) Mapping of variant sites on a kallikrein structure. The overall structure is depicted as a green ribbon. Variant sites are shown as sticks. The catalytic triad and the second SS6 cysteine are shown as lines. Chapter 3 | Papers 99 Figure 4 – Relative abundance of KLK3 p.S210W variant in seminal plasma. Spectral counts for p.S210W residue in two heterozygous (Het_1, Het_2) and of 33 homozygous (Hz) individuals. Total spectral counts are shown for Het_1 and Het_2 individuals, and the mean of spectral counts are displayed for Het_1+2 and Hz. Chapter 4 Final Discussion Chapter 4 | Final Discussion 108 to alterations in selective pressures. Indeed, such fast response patterns are observed in the case of artificial selection in maize (Innan and Kim 2004; Durand et al. 2010). A selective hypothesis associated with reproductive functions was considered at fist given that KLK4 has been also proposed to have a role in the proteolytic cascade of semen liquefaction and considering the reported importance of KLKs in the evolution of semen coagulation in primates (Paper I) (Takayama et al. 2001b). However, male physiological traits indicate that sperm competition is unlikely to have played a part in human evolution and these are rather consistent with a history of monogamy (unimale mating systems) (Dixson 2009). For example, in comparison with other mammals, humans have small testis relative to body size, longer spermatogenesis times, lower sperm cell counts and tinier spermatozoa midpieces. Furthermore, it was hypothesized that the reduced testis size observed in some modern Asian populations is not due to sperm competition, neither to changes in their mating system, but instead, are a byproduct of the selective forces acting on gonadal function in women to prevent multiple ovulation (Dixson 2009). Therefore, taking into account the traits known to have been targeted by natural selection in human populations and the possible implications of KLK4 in these processes, two selective hypothesis related to tooth and epidermal features were advanced (Paper II). On one hand, the expression of KLK4 in some epidermal layers together with its ability to activate in vitro several molecules involved in important pathways of skin physiology, such as keratinization, melanosome transfer and skin desquamation, turned out to be an extremely attractive hypothesis (Seiberg et al. 2000; Komatsu et al. 2003; Babiarz-Magee et al. 2004; Komatsu et al. 2005; Matsumura et al. 2005; Becker-Pauly et al. 2007; Ramsay et al. 2008). Although, protein expression in skin is still a controversial issue and the KLK4 inactivation in amelogenesis imperfecta disease and mouse models does not seem to have major outcomes in epidermal phenotypes (Hart et al. 2004; Simmer et al. 2009; Wang et al. 2013). Nonetheless, it is important to note that skin as well as many human morphological traits, is thought to be polygenic, in which KLK4 may only contribute through a moderate effect and, hence, subtle changes in phenotypes might have remained undetected. On the other hand, KLK4 play an essential role in tooth maturation and accordingly, its disrupted activity result in enamel defects in humans and mice (Hart et al. 2004; Simmer et al. 2009; Wang et al. 2013). Additionally, KLK4 arose by a duplication of KLK5 near the divergence of Boreoeutheria (e.g. primates, rodents and artiodactyls) and Afrotheria (e.g. elephant, hyrax and tenrec) lineages, where the latter taxonomic group presents a characteristic delayed in dental eruption (Kawasaki et al. 2014). Moreover, Chapter 4 | Final Discussion 109 KLK4 has been pseudogenized in toothless minke and bowhead whales (Balaenoptera acutorostrata and Balaena mysticetus, respectively), whereas in the toothed cetaceans (orca - Orcinus orca and bottlenose dolphin - Tursiops truncates) the gene remained intact (Kawasaki et al. 2014; Keane et al. 2015). In humans, and in East Asians in particular, it is possible to observe common dental variations, like upper central incisor shoveling, enamel extensions of the first maxillary molar and other dental traits, often designated as sinodonty, which can be partially associated with the genetic variation in EDAR, another gene showing footprints of positive selection in East Asians (Turner 1990; Hanihara and Ishida 2005; Sabeti et al. 2007; Hanihara 2008; Kimura et al. 2009; Kamberov et al. 2013). For the reasons mentioned above, and bearing in mind that tooth morphologies are also likely a polygenic trait, KLK4 might represent yet another gene contributing to an adaptive tooth phenotype through a moderate effect. Still, the pervasive nature of KLK4, its wide pattern of tissue expression and the observed outcomes of selected variants in three different cellular systems, seem to suggest the existence of pleiotropic effects in other physiological functions, which may result in increased resistance or susceptibility to human disease. 2. Implication of KLKs genetic variation in human health and disease Most, if not all, biological processes in humans comprise complex proteolytic networks, in which spatial and temporal alterations can lead to different disease states (Lopez-Otin and Bond 2008; Quesada et al. 2009). In this scope, the KLK gene family is no exception and despite the numerous studies performed already to understand the impact of KLKs in both normal and pathological conditions, significant pieces of information are still missing. For instance, from a clinical perspective, most of the literature has been focused on cancer risk, where KLK3 is by far the best studied member of the family, due to its early recognition as a biomarker for prostate cancer. Nevertheless, the recent development of diverse Klk mouse models (knockout, knock-in and transgenic) have shed light into the associations of each KLK with distinctive human disorders, such as amelogenesis imperfecta, atopic dermatitis, myocardial ischaemia, multiple sclerosis and schizophrenia (Ny and Egelrud 2004; Pons et al. 2008; Simmer et al. 2009; Smith et al. 2011; Tamura et al. 2012; Murakami et al. 2013; Furio et al. 2014). Nonetheless, the potential implications of KLK variation outside of malignancy remain largely unknown. Chapter 4 | Final Discussion 110 Hence, in-depth comprehensive surveys of the KLK locus within varied diseases, including male infertility, are fundamental. In a third work, the contribution of KLKs variability was explored in the framework of male infertility and three non-mutually exclusive phenotypes, hyperviscosity, asthenozoospermia and oligozoospermia, simultaneously with the study of their substrates and inhibitors, the SEMGs, EPPIN and EPPIN-like genes, all belonging to the WFDC family located at chromosome 20q13 (Paper III). Male infertility is considered a multifactorial disease with a strong genetic component, in which common variants may account with small increments into disease susceptibility (Carrell and Aston 2011). In this particular disorder, it is rather straightforward that strong deleterious alleles will be rapidly wipe out from extant human populations because of their negative impact in the reproductive fitness. Conversely, as mentioned earlier, as human populations have expanded many novel variants were originated and, due to their recent origins, they will still be segregating at low-frequencies in populations in spite of a possible influence in human fertility (codominant inheritance). Altogether, male infertility, or more precisely male subfertility, like other human complex disorders, can be considered as a consequence of a wide range of genetic variants with variable size effects into the different disease phenotypes (Aston and Carrell 2009; Aston et al. 2010; Carrell and Aston 2011; Lettre 2014). Specifically, in the case-control study of male infertility centered on the KLK and WFDC loci, the majority of identified candidate variants with potential serious effects in protein structure, activity and interactions, splicing processing and expression regulation were low-frequency variants, among which KLKs were overrepresented, suggesting a greater burden of these genes in male reproduction in comparison to WFDCs (Paper III). Among the recognized infertility phenotypes, delayed liquefaction, hyperviscosity and asthenozoospermia are the ones more probably correlated to a deregulated activity of KLKs and WFDCs in the semen, given that a small modification of the fine-tuned protein interactions in the semen might potentially interfere in its quality. Consistently, most KLKs have been previously described to be downregulated in individuals exhibiting the hyperviscosity phenotype (KLK1-2, KLK5-8, KLK10, KLK13-14) and a few KLKs were also found to display a significantly reduced protein expression in cases with delayed liquefaction (KLK2-3 and KLK13-14) (Emami et al. 2009). On the other hand, associations between abnormal sperm motility (asthenozoospermia) and KLK14 downregulation, SEMG1 augmented expression and the higher prevalence of EPPIN low-frequency variants were also reported (Martinez-Heredia et al. 2008; Emami et al. 2009; Ding et al. 2010a; Ding et al. 2010b). In overview, the KLK and WFDC genetic variants identified in Chapter 4 | Final Discussion 111 Paper III comprised several SNVs possibly influencing the protein contents in the semen by mechanisms of abnormal protein and mRNA degradation (nonsynonymous substitutions and splice variants) or by gene misregulation (variants located in regulatory regions), which could in some instances relate to the reduced KLK levels already reported. Conversely, several other variants (nonsynonymous substitutions only) were found to likely affect important protein interactions in semen because of the insertion or removal of proteolysis cleavage sites, or even the substrate allowance in the protease catalytic pocket, all of them possibly connected to hyperviscosity and asthenozoospermia phenotypes. Furthermore, the distribution of these susceptibility markers across KLKs seems to suggest KLK3, KLK12 and KLK14 as the genes with greater impact in male fertility. Although a greater focus was given to nonsynonymous replacements and splice variants in the genotyping phases, mainly due to the reduced reliability of bioinformatics predictions for nucleotide substitutions located in transcription factor binding sites, the variants placed at coding regions may only explain a minor part of the genetic predisposition to male infertility (Paper III). Indeed, only a small fraction of the 456 identified variants were nonsynonymous substitutions (12.7%; 9.6% low-frequency), which reinforces a possible contribution of gene regulation in the disease pathogenesis. This finding concurs with previous assumptions made for targets of natural selection, in which eQTLs are also thought to play a significant role in complex diseases by altering the gene expression timings and intensities (Epstein 2009; Nicolae et al. 2010; Bryois et al. 2014). In addition, it has been shown that common SNVs associated with complex traits are more likely to be eQTLs than nonsynonymous variants, in which the effect size is inversely correlated to allele frequency, in brief, the greater the allele frequency, the lower the contribution to the trait (Nicolae et al. 2010; Battle et al. 2014). Interestingly, for a common variant located KLK7 and a low-frequency CNV of SEMG1, both showing significant associations with male infertility, the ancestral alleles were the ones conferring an increased susceptibility to the disease. From an evolutionary perspective, this evidence fits well the “thrifty genotype” theory, which proposes alleles that were once beneficial in past human history, nowadays are less adapted to the environmental conditions of western societies (Neel 1962; Neel et al. 1998). Consistently, other complex disorders linked to this theory, like obesity and diabetes, are known to affect human fertility, as well as lifestyles habits such as smoking and alcohol consumption and tendency for parents having the first child at increased ages (Vine et al. 1994; Sharpe and Franks 2002). Chapter 4 | Final Discussion 112 Finally, the semen is a complex body fluid known to be enriched in proteases and these have been established to have pervasive links with each other, their inhibitors and substrates (Pilch and Mann 2006; Batruch et al. 2011; Laflamme and Wolfner 2013; Fortelny et al. 2014). Therefore, it is attractive to speculate that some of these enzymes may have critical roles in the semen liquefaction profiling with possible functional repercussions in the spermatozoa abilities to fertilize the egg. Globally, these KLK and WFDC variants are far from explaining the genetic basis of male infertility, or even hyperviscosity or asthenozoospermia phenotypes, and only represent a minor proportion of its susceptibility and deleterious variation. In the near future, with the continuous decline of whole exome and genome sequencing costs, together with the advances in proteome technologies, new insights concerning the pathogenesis of male infertility will probably emerge with benefits in patient personalized medicine. Chapter 5 Concluding Remarks Chapter 5 | Concluding Remarks 115 This work provides further support for an implication of proteolysis genes, specifically KLKs, as targets of adaptive evolution in primates through diverse biological functions including male reproduction and for a possible role of KLK sequence variation in both beneficial traits and disease phenotypes. 1. It sustains the hypothesis that the evolution of KLK2 and KLK3 might have been driven by reproductive biology in a sperm competition manner. This study contributed not only to consolidate the origin of KLK3 in Catarrhini through an event of duplication and functional divergence from KLK2 towards a chymotrypsin-like activity and, consequently, an extended enzymatic spectrum over SEMG1 and SEMG2 cleavage, but also to unveil a complex dynamics of KLK2 and KLK3 evolution mediated by different genomic mechanisms of gene loss in a close correlation to SEMG gene structure, primate mating system and semen coagulation rates. 2. It clarifies a signal of natural selection shaping the genetic diversity of the KLK cluster in Asian populations. A complex signature of positive selection was disclosed in East Asians favoring a haplotype defined by three variants (rs1654456_G, rs198968_T and rs17800874_A) acting synergistically to downregulate KLK4. A reasonable hypothesis is that this haplotype was driven to high frequencies in East Asians by offering a selective advantage into phenotypic traits characteristic of these populations, such as tooth shape and structure, and epidermal features. Still, due to the pervasive nature of KLK4, the possibility of pleiotropic effects in other physiological functions, including in reproductive biology, also exists and may result in an increased disease resistance or susceptibility. In a wider perspective, this study endorses the concept that soft sweeps and polygenic adaptation may actually be more common in human evolution than classic selective sweeps. 3. It highlights several low-frequency variants at KLK and SEMG genes that may contribute to different male infertility phenotypes. An excess of low-frequency variants was observed for the KLK cluster but not for the WFDC locus and this burden was independent of the infertility phenotype considered. Among the most promising variants, two were found to overlap onto the protein structure in KLK3 (p.E131K) and in KLK14 (rs112658494, p.R138W). Furthermore, one variant in KLK12 (rs61742847, p.P34L) and another in SEMG1 (rs147894843, p.G400D) were significantly associated with semen hyperviscosity and asthenozoospermia, respectively. In addition, a common Chapter 5 | Concluding Remarks 116 intronic variant in KLK7 (rs1654526) and the copy number variation embracing one of the repeat units in SEMG1 were connected to a protective role against hyperviscosity and oligozoospermia, respectively. Nevertheless, the identified candidate variants are expected to only play a part in the genetic background of male infertility, given that they were always found in heterozygosity and that KLKs often overlap in substrate affinity. Chapter 6 References