On the tear proteome of the house mouse (Mus musculus musculus) in relation to chemical signalling

Ethical standards

All animal procedures were carried out in strict accordance with the law of the Czech Republic paragraph 17 no. 246/1992 and the local ethics committee of the Faculty of Science, Charles University in Prague chaired by Dr. Stanislav Vybíral specifically approved this study in accordance with accreditation no. 27335/2013-17214 valid until 2019. Animals were sacrificed by cervical dislocation.

Animals

Fourteen individuals of the house mouse (the eastern form, M. m. musculus) used in this study were captured in the Czech Republic near Bruntál - 49.9884447N, 17.4647019E (one male, one female), in Velké Bílovice - 48.8492886N, 16.8922736E (three males, three females), Prague-Bohnice - 50.1341539N, 14.4142189E (three males, three females). All animals were trapped in human houses and garden shelters. On the day of capture or the next day, all animals were transferred to our animal facility. Each animal was caged individually with ad libitum access to water and food.

Tear collection

Eye lavage was used as a non-invasive method of tear collection. Both eyes were carefully rinsed with 10 µl of the saline physiology solution by a gentle pipetting and samples were then pooled. The process was repeated three times with at least a two hour interval between every rinsing, and each sample was analysed twice with MS to produce mean values from the methodology duplicates. This was done in the ‘in-house’ Mass Spectrometry and Proteomics Service Laboratory, Faculty of Science, Charles University in Prague.

Protein digestion

Protein samples were precipitated with the ice-cold acetone and followed by a re-suspension of dried pellets in the digestion buffer (1% SDC, 100 mM TEAB—pH  = 8.5). Protein concentration of each lysate was determined using the BCA assay kit (Fisher Scientific, Waltham, MA, USA). Cysteines in 20 µg of proteins were reduced with a final concentration of 5 mM TCEP (60 °C for 60 min) and blocked with 10 mM MMTS (i.e., S-methyl methanethiosulfonate, 10 min Room Temperature). Samples were cleaved with trypsin (i.e., 1/50, trypsin/protein) in 37 °C overnight. Peptides were desalted on a Michrom C18 column.

nLC-MS² analysis

Nano Reversed phase columns were used (EASY-Spray column, 50 cm × 75 µm ID, PepMap C18, 2 µm particles, 100 Å pore size). Mobile phase buffer A was composed of water, 2% acetonitrile and 0.1% formic acid. Mobile phase B contained 80% acetonitrile, and 0.1% formic acid. Samples were loaded onto a trap column (Acclaim PepMap300, C18, 5 µm, 300 Å Wide Pore, 300 µm × 5 mm, five Cartridges) for 4 min at 15 µl/min loading buffer was composed of water, 2% acetonitrile and 0.1% trifluoroacetic acid. After 4 min ventile was switched and Mobile phase B increased from 2% to 40% B at 60 min, 90% B at 61 min, hold for 8 min, and 2% B at 70 min, hold for 15 min until the end of run.

Eluting peptide cations were converted to gas-phase ions by electrospray ionization and analysed on a Thermo Orbitrap Fusion (Q-OT-qIT; Thermo Fisher, Waltham, MA, USA). Survey scans of peptide precursors from 400 to 1,600 m/z were performed at 120K resolution (at 200 m/z) with a 5 × 10⁵ ion count target. Tandem MS was performed by isolation at 1.5 Th with the quadrupole, HCD fragmentation with normalized collision energy of 30, and rapid scan MS analysis in the ion trap. The MS² ion count target was set to 10⁴ and the max injection time was 35 ms. Only those precursors with charge state 2–6 were sampled for MS². The dynamic exclusion duration was set to 45 swith a 10 ppm tolerance around the selected precursor and its isotopes. Monoisotopic precursor selection was turned on. The instrument was run in top speed mode with 2 s cycles.

Protein analysis

All data were analysed and quantified with MaxQuant software (version 1.5.3.8) (Cox et al., 2014). The false discovery rate (FDR) was set to 1% for both proteins and peptides and we specified a minimum peptide length of seven amino acids. The Andromeda search engine was used for the MS/MS spectra search against the Uniprot Mus musculus database (downloaded on June, 2015), containing 44,900 entries. Enzyme specificity was set as C-terminal to Arg and Lys, also allowing cleavage at proline bonds (Rodriguez et al., 2008) and a maximum of two missed cleavages. Dithiomethylation of cysteine was selected as fixed modification and N-terminal protein acetylation and methionine oxidation as variable modifications. The “match between runs” feature of MaxQuant was used to transfer identifications to other LC-MS/MS runs based on their masses and retention time (maximum deviation 0.7 min) and this was also used in all quantification experiments. Quantifications were performed with the label-free algorithms described recently (Cox et al., 2014) using a combination of unique and razor peptides. To detect differentially expressed / abundant proteins, we used the Power Law Global Error Model (PLGEM) (Pavelka et al., 2004) within the Bioconductor package in R software (Gentleman et al., 2004). This model was first developed to quantify microarray data (Pavelka et al., 2004); however, due to similar statistical properties—namely the distribution of signal values deviating from normality—it has proved to be an amenable model for the quantification of label-free MS-based proteomics data (Pavelka et al., 2008). Next, we calculated the signal-to-noise ratio—STN (equation provided in Pavelka et al., 2008), because it explicitly takes unequal variances into account and because it penalizes proteins that have higher variance in each class more than those proteins that have a high variance in one class and a low variance in another (Pavelka et al., 2004). PLGEM was fitted on a set of replicates from female data, thus setting experimental baseline. Correlation between the mean values and standard deviations was high (r² = 0.96, Pearson = 0.94) so we continued with the resampled STNs and calculated differences with corresponding p-values between males and females.

Protein surface modelling

The surface electrostatics modelling involved several steps. First, we downloaded the structures from the RSCB Protein Data Bank (http://www.rcsb.org/) under accession IDs: 3S26, 1I04 and 2L9C, respectively. Because the mouse OBP1 structure has no record in the database we had to predict it with i-TASSER (Iterative Threading ASSEmbly Refinement) program (http://zhanglab.ccmb.med.umich.edu/I-TASSER/) with the rat ortholog OBP1F (PDB ID: 3FIQ, 76% similarity) as template for the homologous modelling. Next, we used PyMOL - Molecular Graphic System (version 1.7.0.0) with APBS (Adaptive Poisson-Boltzman Solver) plugin to model the electrostatics with the default software settings.

RNAseq: samples, cDNA, sequencing, and analysis

Individual mice were sacrificed next day after the last tear sampling. The exorbital lacrimal glands (i.e., one gland per individual, the other one is stored) were dissected and immediately placed into RLT buffer (Qiagen, Hilden, Germany) and homogenised in MagNALyser (Roche) for 30 s at 6,000 rpm. RNA was isolated using the RNeasy Mini Kit (Qiagen, Hilden, Germany) according to the manufactures protocol with on-column DNase I treatment. The purity and concentration of eluted RNA was measured with a NanoDrop ND1000. The quality of RNA was checked on agarose gel electrophoresis (AGE). RNA was stored at −70 °C pending further use. For the next step, we selected only high quality samples from four male and four female replicates.

cDNA was prepared using the SMARTer PCR cDNA Synthesis Kit (Clontech, Mountain View, CA USA) and amplified with Advantage 2 PCR Kit (Clontech, Mountain View, CA USA). Both procedures were handled according to protocol for Trimmer-2 Normalization Kit (Evrogen, Moscow, Russia). The products of optimalized cDNA amplification were then loaded on AGE. For each sample, only the area of product in range from ∼400 bp to ∼1,300 bp (well visible area full of bands) was excized from the gel and the DNA products were extracted using the Gel/PCR DNA Fragments Extraction Kit (Geneaid, New Taipei City, Taiwan). Appropriate amounts of size-selected products were then secondarily amplified according to the recommended protocol from Evrogen. Products of secondary amplification were purified using MiniElute PCR Purification Kit (Qiagen, Hilden, Germany). Purified products (and the range where they emerge) were checked on AGE. Purity was analysed with NanoDrop ND1000. Concentration was measured/determined using Quant-it Pico Green dsDNA Assay Kit (Invitrogen, Carlsbad, CA, USA) and fluorimeter (Hoefer DQ 300). Rapid Library (RL) was prepared for each transcriptome (four males and four females) according to Rapid Library Preparation Manual (my454.com). Equal amounts from each of eight Rapid Libraries (10⁷ molecules per µl dilution) were mixed and then used for emPCR.

We conducted 454 RNA-sequencing with a desktop pyro-sequencer GS Junior from Roche using the long reads mode. To increase the precision of transcript mapping, we excised from a gel and sequenced only transcripts between ∼400 and 1,300 bp. Transcripts of this length include those of genes, described for their involvement in chemical communication (e.g., lipocalins). This method is amenable to further analyses because the nebulization step is skipped and, therefore, whole transcripts instead of their fragments are further pyro-sequenced and mapped. We estimated particular expression levels from the number of uniquely mapped transcripts assigned to each annotated gene. All steps followed the provider’s instructions for sequencing with GS Junior (emPCR Amplification Method Manual Lib-L and Sequencing Method Manual; Roche, Basel, Switzerland). We obtained >165,000 high quality (HQ) reads. HQ 454 Reads were multiplexed, trimmed (i.e., using a trimming database that contains primers used for library preparations), filtered and aligned into contigs against Mus musculus cDNA database (“the super-set of all known, novel and pseudo gene predictions”; ensembl.org, 17-FEB-2015 version) and using GS Reference Mapper (Roche, Basel, Switzerland). Differential expression was analysed in R software using the DEseq routine within the Bioconductor package (Gentleman et al., 2004).

RNA-seq data availability

The transcriptome data is provided as bam files in ‘Sequencing Read Archive’ (www.ncbi.nlm.nih.gov/sra) under the accession numbers SRP063762 and BioProject PRJNA295909.

The tear proteome and the level of sexual dimorphism

We have generated the tear proteome of the house mouse, M. m. musculus and detected a total of 719 proteins at 0.01 FDR (i.e., False Discovery Rate for all peptides and proteins). First of all, we reduced our data such that only the proteins that were detected in three or more individuals were further analyzed (i.e., 457 proteins). Our aim was to identify those proteins that are sexually dimorphic (Fig. 1) and those that represent the top 5% of the most abundant proteins that may characterize the mouse tears (Fig. 2A).

PLGEM analysis of the level of sexual dimorphism revealed that 68 (14.9%) out of 457 proteins identified at 1% FDR and p < 0.05 were sexually dimorphic, Fig. 1. Male biased proteins included 36 (7.8%) and female biased proteins included 32 (7%) successful identifications (i.e., listed in Data S1). Thus, male-biased proteins were not more common than female-biased proteins in the tear proteome. The most significant dimorphic proteins (i.e., top 5% in Fig. 2B) included the female-biased OBP7, the male-biased MUP4, the male-only ESP1, the male-biased ESP38, and several male-biased secretoglobins (SCGB1B19, SCGB1B3, SCGB2A2—Mammaglobin, SCGB2B3, SCGB2B7). Secretoglobins are found in mammalian secretions and have important roles in the modulation of inflammation and tissue repair (Jackson et al., 2011), are involved in removing toxins (Zhou et al., 2011), and may have roles in chemical communication (Karn & Laukaitis, 2015). Kallikrein 1-related peptidases were also significantly sexual dimorphic (i.e., female-biased KLK1B22, KLK1B1, and KLK1B3), however, this pattern (though significant) was not consistent across all the females tested. KLKs are known as network members that are crucial for homeostasis of stratified epithelia and for activating antimicrobial Cathelicidins (Kasparek et al., 2017). Interestingly, we have also detected sexually dimorphic BPI proteins. Bactericidal/permeability-increasing proteins (BPI) are ∼50 kDa proteins that are a part of the innate immune system, and have an antibacterial activity against the gram-negative bacteria (LeClair, 2003b). We have detected three BPIs, of which BPIFA2 was male biased, BPIFA6 was female biased, whilst males and females equally expressed BPIFB9B.

The most abundant tear proteins

Based on the median value we sorted our data to detect the most abundant proteins in the tear proteome. The top 5% of the most abundant proteins that characterize the soluble tear proteome of the mouse are depicted in Fig. 2A, and include for example the female-biased lipocalins OBP5, OBP7, the unbiased lipocalins OBP1 and LCN11, and the male-biased lipocalin MUP4, Fig. 2C. Other proteins dominating the soluble tear-proteome included three male-biased secretoglobins (SCGB1B3, SCGB1B20/SCGB1B29, SCGB2B20/SCGB2B7), two unbiased secretoglobins (SCGB1B2, SCGB2B2), male-biased carbonic anhydrase 6 (CAH6), (unbiased) exocrine secreted peptide ESP6, Lacrein, and female-biased prolactin inducible protein (PIP). Interestingly, out of the top 5% most abundant proteins, ∼50% of them (i.e., 12) were significantly sexual dimorphic. Thus, even though the level of sexual dimorphism is rather low within the complete tear proteome (i.e., 15%), those few proteins that were most abundant were often the most sexually dimorphic.

Sex-unique proteins

We provide visual representation of all proteins using MA plot, also including potentially sex-unique proteins (Fig. 1), where significant points are colored from green (p < 0.05) to blue (p < 0.01). Female-unique proteins included S10A8/S10A9 which are calcium- and zinc-binding proteins and which play important roles in the regulation of inflammatory processes and immune responses, and can induce neutrophil chemotaxis and adhesion (Vogl et al., 2007). We have also detected the female-unique secretoglobin SCGB2B20. Other female-unique proteins invlolved RENI2, LIPR1 and one keratin (KT33A). Male-unique proteins included the Secretoglobin SCGB1B19, the exocrine gland-secreted peptide ESP1, Zn-Alpha2-Glycoprotein (i.e., ZA2G), zona pellucida-like domain-containing protein 1 (CUZD1), and products of the two predicted genes Gm12887 and Gm1330. ESP1 was co-expressed with other ESPs (ESP3, ESP4, ESP6, ESP15, ESP16, ESP18, ESP34, ESP38) in tears, Fig. 2D. Though male-unique in tears, ESP1 is present in male and female saliva (Stopka et al., 2016). To add, various proteins that were previously detected as sex-unique now seem to be rather sex-biased and not sex-unique/limited when new LC-MS/MS techniques with higher detection limits are employed instead of the gel-based MS techniques (e.g., Karn & Laukaitis, 2015).

Transcriptome: mRNAseq based analysis of exorbital lacrimal glands

To detect lipocalins that are excreted from the exorbital lacrimal glands, we performed the analysis of size-selected transcriptomes from ∼400 bp to ∼1,300 bp demonstrated with histograms in Fig. 3A. Top ten percent of highly expressed transcripts included Scgb2a2, Pip, Lcrn, Scgb2b2, Scgb1b2, Scgb2b24, Spt1, Lcn11, Obp5, Scgb2b7, Esp6, Esp15, Scgb1b3, Sval2, Obp7, Scgb2b3, Bpifa2, Bglap3, Esp16, Scgb1b7, Gm20594 (Mtrnr2l), Wfdc18 (underlined are the transcripts encoding proteins that were detected within the top 10% of highly expressed proteins). However, direct comparison between transcriptomic and proteomic datasets was not possible, because the transcriptome was prepared from selected ranges of gel extracted tissue mRNAs whilst the tear proteome over-represents the soluble fractions of all proteins. Next, we searched for sexually dimorphic genes that may account for sex-specific differences with the DESeq routine within the Bioconductor package (Gentleman et al., 2004). We have filtered for further analysis only the data where the sum of counts per row ≥10. Then, we normalised the data with a size factor vector to make the libraries comparable. Because DESeq calculates sexual dimorphisms from the original non-transformed number of counts we first looked at the level of variation between replicates within sex. When dispersion values are plotted against the means of the normalised counts (Fig. 3B) it is evident from the slope of the red fitting curve that data with a low mean of normalized counts have higher levels of dispersion than high expression data.

Next, we searched for differentially expressed genes by calling the nbinomTest in DESeq, vizualized in Fig. 3C. We have detected a total of 6 female-biased genes (Obp5, Obp7, Obp8, Spt1, Hba, and Scgb2b1) and a total of 17 male-biased genes (Scgb2b7, Scgb1b20, Scgb1b3, Scgb1b7, Esp18, 9530002B09Rik, Scgb2b3, Esp16, Mup4, Esp24, Cyp4a12b, RP23-421B1.4, Nop10, Scgb1b28-ps, Esp1, Scgb2b20, and Pigr), Fig. 3D. Next we asked which of the above sex-biased genes are most differentially expressed. Using Benjamini–Hochberg corrections we have generated new ‘p-adjusted’ values. These genes (p-adjusted < 0.05) included a total of 13 genes with female-biased Obp5, Obp7, Obp8, and Spt1, whilst male-biased genes included Mup4, five secretoglobins (Scgb2b7, Scgb1b20, Scgb1b3, Scgb1b7, and Scgb2b3), two Esp s (Esp16, Esp18) and the gene 9530002B09Rik (synonym: Vpp1—Ventral prostate predominant l, which was originally thought to be exclusively expressed in the prostate (Wubah et al., 2002)). The resulting pattern is plotted using MA plot (Fig. 3C) with red colouring of those genes that are significant at p-adjusted < 0.05 whilst all data with p < 0.05 (not corrected) are vizualized in Fig. 3D and those that have p-adjusted <0.05 are depicted with asterisks. Partial support for mRNA/protein concordance is provided in Fig. 3E, showing that those transcripts that are significantly sexually dimorphic and represent the soluble protein fraction are also detected on the level of protein. Furthermore, significant sexual dimorphisms of Obp7 and Mup4 (Figs. 3D, 3E) are in agreement with our proteomic analysis in this study and with our previously published study using qPCR (Stopkova et al., 2016). Our simple RNAseq analysis revealed the expression of Obp8 which was previously detected only with bioinformatics tools (Stopková et al., 2009; Stopkova et al., 2016), thus providing the first evidence for the expression of Obp8 transcript. Interestingly, OBP1 was one of the most abundant tear proteins in individuals of both sex. However, the expression of Obp1 transcript was medium/low in this study and low with qPCR in our previous study (Stopkova et al., 2016). Thus, OBP1 could be a product of several other tear-secreting glands including infra-orbital glands, accessory lacrimal glands and/or epithelial cells of ocular mucosa.

Anti-microbial peptides

BPI proteins have an antibacterial activity against gram-negative bacteria (LeClair, 2003b). The saliva proteome contains seven members of the bactericidal/permeability-increasing proteins (i.e., BPI Leclair, 2003a; LeClair, 2003b) which are male biased (Stopka et al., 2016) and include BPIA1, BPIB1, BPIB2, BPIB3, BPIFA2, BPIFB5, BPIFB9B (Stopka et al., 2016). However, tears only contain BPIFA2/Bpifa2—one of the most expressed transcripts, female-biased BPIFA6, and un-biased BPIFB9B. Thus, we searched for other proteins/peptides which may have similar roles due to their amphipathic structural properties or proteolytic activities. Recently, WFDC proteins (i.e., ‘Whey acidic proteins four disulphide core’) were shown to have anti-microbial properties (Scott, Weldon & Taggart, 2011) and the two members WFDC12 and WFDC18 are present in mouse saliva as proteins encoded by submandibular gland transcripts (i.e., Wfdc12, Wfdc18) (Stopka et al., 2016). In this study, we have detected WFDC12 and WFDC18 as transcripts of the exorbital lacrimal glands (i.e., Wfdc12, and the highly expressed Wfdc18), but only WFDC18 was detected in tears on the proteomic level and just in two males. Our results, however, provide evidence that the major antimicrobial protein in tears is TRFL (Lactotransferrin), Fig. 2A. Lactotransferrin also known as lactoferrin (LF) has antimicrobial properties (bactericidal, fungicidal) and is a part of the innate immune system, mainly at mucoses (Sanchez, Calvo & Brock, 1992). In the tear proteome, we detected the male-biased TRFL as one of the most abundant proteins and similar amounts were previously also detected in saliva (Stopka et al., 2016).

Homology modelling of representative lipocalin structures

OBPs and MUPs are likely to have complementary roles because OBPs are less hydrophobic and have higher iso-electric points whilst MUPs are more acidic and hydrophobic (Stopkova et al., 2016). In Fig. 4, we provide four representative lipocalin structures from homology modelling, thus showing that different lipocalins have similar structures but different electrostatics properties. The distribution of negative and positive residua is not random in OBP1 and even less so in LCN2. The structure of LCN2 is amphipathic because the interaction is driven by ionic strength between positively charged amino acid residues at the barrel opening of LCN2 and negatively charged siderophores. Thus, LCN2 is antimicrobial, as it efficiently scavenges for catecholate-type siderophores which bacteria produce to scavenge for free iron (Flo et al., 2004). The structure of OBP1, however, is amphipathic due to a non-random distribution of positively charged residua on its surface and may potentially be antimicrobial (i.e., similar to CRAMP/CAMP (Gallo et al., 1997)). In addition, such amphipathic structure of OBP1 may also aid to a direct attack upon negatively charged bacterial membranes by its oppositely charged OBP1 surface residua including the positively charged alpha helix.