Age- and Sex-Adjusted Reference Intervals in Tear Cytokine Levels in Healthy Subjects

Fernández, Itziar; Enríquez-de-Salamanca, Amalia; Portero, Alejandro; García-Vázquez, Carmen; Calonge, Margarita; Herreras, José M.

doi:10.3390/app11198958

Open AccessArticle

Age- and Sex-Adjusted Reference Intervals in Tear Cytokine Levels in Healthy Subjects

by

Itziar Fernández

^1,2,3

,

Amalia Enríquez-de-Salamanca

^1,3,*

,

Alejandro Portero

^1,†,

Carmen García-Vázquez

¹

,

Margarita Calonge

^1,3

and

José M. Herreras

^1,3,4

¹

IOBA (Institute of Applied Ophthalmobiology), University of Valladolid, E-47011 Valladolid, Spain

²

Department of Statistics and Operations Research, Universidad de Valladolid, E-47011 Valladolid, Spain

³

Biomedical Research Networking Center in Bioengineering, Biomaterials and Nanomedicine (CIBER-BBN), E-47011 Valladolid, Spain

⁴

Ophthalmology Service, University Clinic Hospital, E-47003 Valladolid, Spain

^*

Author to whom correspondence should be addressed.

^†

Present address: Immunology and Ocular Surface Unit, Hospital la Zarzuela, E-28023 Madrid, Spain.

Appl. Sci. 2021, 11(19), 8958; https://0-doi-org.brum.beds.ac.uk/10.3390/app11198958

Submission received: 25 June 2021 / Revised: 17 September 2021 / Accepted: 23 September 2021 / Published: 26 September 2021

(This article belongs to the Special Issue Potential Biomarkers in Tears)

Download

Browse Figures

Versions Notes

Abstract

:

Featured Application

Age- and sex-adjusted reference intervals (RIs) could be an important tool to extend the utility of tear cytokine levels as biomarkers of ocular inflammatory diseases. Particularly, age- and sex-adjusted RIs established in this study for tear levels in healthy subjects have shown effectiveness for severe DED diagnosis.

Abstract

Alterations in tear cytokine levels have been associated with various ocular disorders as compared to those in healthy subjects. However, age and sex are not always considered in these comparisons. In this study we aimed to establish age and sex reference intervals (RIs) for tear cytokine levels in healthy people. Tear samples were taken from 75 males and 82 females, aged 18–88 years, and tear cytokine levels were determined. Age- and sex-adjusted RIs for epidermal growth factor (EGF), fractalkine, interleukin (IL)-1 receptor antagonist (RA), IL-7, IL-8, interferon inducible protein (IP)-10, monocyte chemotactic protein (MCP)-1, and vascular endothelial growth factor (VEGF) tear cytokine levels in a healthy sample were established using generalized additive for location, scale and shape (GAMLSS) models. RIs were tested in two external samples: a validation sample of 40 individuals with normal results at four Dry Eye Disease (DED) clinical diagnostic tests (OSDI, T-BUT, corneal staining and Schirmer test); and a utility sample of 13 severe DED cases. IL-1RA, IL-8, IP-10, and MCP-1 levels showed a positive association with age, while EGF was negatively correlated. IL-7 concentration increased up to 40 years and again after 70 years, observing a quasi-linear decrease between them. For VEGF, higher levels were observed in the middle-aged range. Regarding sex-influence, fractalkine tear levels were higher in men, whereas those of IL-7, IL-8, and IP-10 were higher in women. Using the estimated age- and sex-adjusted RIs, more than 92% of the validation sample was correctly classified, and 100% of the severe DED patients in the utility sample had concentrations outside the RIs in at least two of the cytokines evaluated.

Keywords:

reference intervals; tear cytokines; age- and sex-influence; GAMLSS models; dry eye disease; biomarkers

1. Introduction

Cytokine profiles have been proposed as biomarkers for numerous infectious and chronic diseases due to their pro- and anti-inflammatory effects [1]. Numerous studies have also reported important alterations in tear cytokines associated with various ocular disorders, such as primary open-angle glaucoma, ocular chronic graft versus host disease, uveitis, ocular allergies, and dry eye disease (DED), among others [2,3,4]. In particular, a highlight is the usefulness of these molecules as biomarkers for DED [5]. For instance, based on cytokine profiles, diagnostic, disease, severity, and therapeutic biomarkers have been established for DED by several authors [6,7,8]. It is widely accepted that ocular surface inflammation has a central role in the occurrence and severity of the DED condition [9,10]. Therefore, a useful tool in evaluating the inflammatory profile of DED may be to analyze the inflammatory mediators in tears, such as tear cytokine profile. In addition, unlike other biological fluids used for cytokine profiling (such as pleural fluid, serum, or cerebrospinal fluid) tear collection is considerably less expensive and requires minimally invasive procedures. In addition, for ocular surface diseases, tear fluid is close to the disease site so it would be expected that tear cytokine profile would be an effective biomarker of them [11].

However, despite its advantages, the use of tear cytokine profiles in clinical practice is not yet widely used as a diagnostic tool. The most common way to evaluate the potential diagnostic value of cytokines has been to establish differences between two or more groups of patients and/or controls. These studies investigate the ability of cytokines to distinguish between the evaluated conditions and confirm the diagnosis. However, a more general approach would be to build reference intervals (RIs) that allow establishing “normal” versus “abnormal” cytokine levels from a sample of healthy subjects. RIs provide a quantitative diagnostic rule to classify a value as potentially pathological when this value could be considered as extreme. The value will be considered extreme if there is only a small proportion of subjects in the healthy population with larger or smaller values. Therefore, a critical issue will be collection of a representative sample from a healthy population. Unfortunately, the number of studies specifically dedicated to investigate tear cytokine profiles in healthy subjects is limited [12,13,14].

Commonly, a RI represents the interval between a pair of symmetrically placed extreme percentiles. Both the 95% RIs (2.5th–97.5th percentile) and the 90% range (5.0th–95.0th percentile) are in use in the literature. Consequently, the RI will be a global range and any patient outcome outside this range will be considered “abnormal”. However, global reference values will be no useful since demographic characteristics, such as age and sex, have been shown to have an influence on biological variations of some protein and cytokine levels in normal human tears [11,12,15,16]. In this research we used generalized additive models for location, scale and shape (GAMLSS) models in a large sample size of subjects without ocular pathology in order to know if age and sex are determining factors and to build age- and sex-adjusted RIs for several cytokine concentrations in tear samples. GAMLSS are very flexible semi-parametric regression models that allow modeling all the distributional parameters: location, scale and shape (skewness and kurtosis), as linear, non-linear parametric, or smooth non-parametric functions of explanatory variables [17]. Assumptions on a distribution form for the response variable are made, allowing a variety of different distribution families. This regression technique provides a platform to fit, compare and check many different models. GAMLSS methodology has been adopted by World Health Organization for creating reference growth curves [18], and has proven to be a useful and flexible technique to estimate adjusted RI (see for example [19,20] among others). In addition, to evaluate the degree to which these established RIs were able to diagnose DED two external samples, one of individuals without DED and another of patients with severe DED, were used.

2. Materials and Methods

2.1. Study Sample

This study was approved by the Institute of Applied Ophthalmobiology (IOBA) institutional review board and by the Ethics Committee of the University of Valladolid, and complied with the Tenets of the Declaration of Helsinki.

Healthy adults aged 18–88 years old with no eye diseases fulfilling the inclusion criteria and that signed the informed consent were included in this study. All subjects were evaluated by the same ophthalmologist (AP) after a complete systemic clinical history. Inclusion criteria were: men or women aged ≥18 years without inflammatory signs in the ocular surface at least in the 12 months prior to inclusion. Exclusion criteria were the use in the 3 previous months of (1) the following systemic medications: immunosuppresants, anti-inflammatory treatment, antidepressants, neuroleptics, anti-hypertensives, corticosteroids, hormones, vitamins, diuretics, antihistamines, beta-blockers, anticholinergics, and retinoids and (2) the use of any topical therapy and artificial tears or lubricants; (3) any type of intraocular, scleral, corneal, and/or conjunctival surgery, either with refractive or therapeutic purposes; (4) any inflammatory pathology with ocular involvement (uveitis, diabetes, allergy, systemic diseases with ocular involvement, collagenopathies); (5) contact lens users of any type in the last year.

2.2. Tear Collection and Tear Inflammatory Molecule Analysis

A non-stimulated basal tear fluid sample (4 µL) was collected from the external canthus of the eye using a 4 µL calibrated glass microcapillary tube (Cat. 1-000-0040, Drummond Scientific, Broomall, PA, USA) in a non-traumatic manner, avoiding reflex tearing as much as possible, as previously described [21,22], until the capillary tube was completely full of liquid. Samples were diluted 1/5 in a 0.5 mL microtube (Sarstedt AG&Co, Nümbrecht, Germany) containing ice-cold Cytokine Assay Buffer (Merck Millipore, Millipore Iberica, Madrid, Spain) and frozen at −80 °C until analysis.

Eighteen cytokines in tear samples were quantified using a Luminex IS-100 equipment (Luminex Corporation, Austin, TX, USA) with Milliplex commercial kits HCYTO60K-17-plex and MPXHCYP2-62K-IL-23 single plex, (both from Merck-Millipore, Millipore Iberica, Madrid, Spain). The concentrations of epidermal growth factor (EGF), fractalkine, interferon (IFN)-γ, interleukin (IL)-10, IL-12p70, IL-15, IL-17A, IL-1β, IL-1 receptor antagonist (IL-1RA), IL-2, IL-6, IL-7, IL-8, interferon gamma-induced protein 10 (IP-10), monocyte chemotactic protein (MCP)-1, tumor necrosis factor (TNF)-α, vascular endothelial growth factor (VEGF), and IL-23 were analyzed according to the manufacturer’s instructions. Briefly, 10 μL of the 1/5 diluted sample were incubated under agitation overnight at 4 °C with beads coated with antibodies specific for each molecule. After washing, the beads were incubated with biotinylated human antibodies for 1 h, followed by incubation with streptavidin-phycoerythrin for 30 min. Standard curves were used to convert fluorescence units to concentration units (pg/mL). Data were stored and analyzed with the “Bead View Software” (Upstate-Millipore Corporation, Watford, UK).

Regression on order statistics (ROS) method [23] was used to impute values below the limit of detection. Molecules that were detected in less than 50% of the samples were not further analyzed. Levels of all other cytokines were analyzed as base 2 log-transformed variables.

2.3. Data Analysis

Quantitative variables were expressed as mean ± standard deviation (SD). Percentages and theirs 95% confidence intervals (95% CI) were used to summarize distributions of qualitative variables. Cytokine levels were described by sex and age groups, considering four equally sized groups of age. Data statistical analysis (including ROS method to impute values below the limit of detection; descriptive statistics; fit, compare and check GAMLSS models; and the use of the RI in external samples) was performed using R Statistical Software (Foundation for Statistical Computing, Vienna, Austria).

2.3.1. Reference Intervals Determination

We calculated percentile curves for each one of the considered cytokines as a function of age stratified by sex using generalized additive models for location, scale and shape (GAMLSS) [17]. GAMLSS are very general and flexible methods for modeling a response variable, in our case a cytokine level, as a function of explanatory variables: age and sex. Any parametric distribution, including highly skewed and kurtotic distributions, can be fitted to the response variable by modeling four parameters: a location parameter, denoted by μ; a scale parameter (σ); and two shape parameters: ν and τ modeling skewness and kurtosis, respectively. Each of the parameters can be a parametric or a non-parametric smooth function of the explanatory variables. In this way, GAMLSS allow to overcome some of the generalized linear models’ limitations. First, the distributional assumption for the response variable is relaxed and a general distribution family can be used. In addition, we can model as a function of explanatory variables, not only the location parameter but also other parameters of the distribution.

Different distributions were fitted to the observed distribution of each molecule concentration on log2 scale. Four distribution families were considered: normal; Box-Cox Cole and Green to explore the existence of skewness in the data; power exponential to explore the possibility of kurtosis; and Box-Cox power exponential distribution, to show if both skewness and kurtosis are present in the data. The age influence on parameters of the considered distributions was modeled as a constant, as a linear function or as a cubic spline. The sex was added to the model as qualitative variable. For a particular cytokine, Akaike Information Criterion (AIC) [24] was used in a stepwise algorithm to choose the best model. In the framework of each distribution family, the strategy for selecting relevant explanatory variables for each of the parameters that characterize the distribution is summarized in Figure 1. Finally, the best model was chosen by selecting the distribution whose final model was better fitted to the data.

2.3.2. Model Diagnostics

To check the adequacy of GAMLSS fitted models we used normalized (randomized) quantile residuals [25]. When the model is correct, residuals will have a standard normal distribution, with mean 0, variance 1, and skewness and kurtosis 0 and 3, respectively. Each fitted model was checked by using summary statistics of the quantile residuals and four graphical tools: residuals against the fitted values of the μ parameter; residuals against each one of the included explanatory variables; a kernel density estimate of the residuals; and the QQ-normal plot of the residuals. In addition, worm plots [26] were used as a diagnostic tool to assess if adjustments for some of fitted model parameters were required (see Residual Analysis in Appendix A for complete information).

2.3.3. Reference Intervals Utility

In order to know the utility of the age- and sex-related RIs estimated with GAMLSS model, we used them to classify subjects from two external independent samples, hereinafter referred to as validation and utility samples, obtained from our own database. Individuals with tear cytokine concentration data along with ocular surface clinical parameters related to the DED diagnosis were selected.

The validation sample consisted of 40 healthy individuals with a mean age of 31.74 ± 15.85 years. There were 19 (47.5%) males and 21 (52.5%) females. Their respective age was 39.31 ± 17.92 years for males and 24.90 ± 9.88 years for females. This healthy sample had normal results at four DED diagnostic tests: the ocular surface index (OSDI) score ≤ 12 points, fluorescein tear break-up time (T-BUT) ≥ 7 s, corneal fluorescein staining (Oxford scale) = 0, and Schirmer test without topical anesthesia ≥ 5 mm in 5 min [27]. The subjects of this sample should be classified as “normal” subjects by the RIs.

By contrast, the utility sample was used to evaluate the capability of these RIs for the diagnosis of DED. Consequently, tear molecule data from 13 cases of severe DED were selected. There were 4 (30.8%) males and 9 (69.2%) females. Their age was 58 ± 7.39 years for males and 58.33 ± 9.72 years for females. The clinical selection criteria for these subjects were having an OSDI score ≥ 33 points, T-BUT < 3 s, corneal fluorescein staining ≥ 3, and Schirmer test < 3 mm in 5 min [27]. Tear sample collection and tear cytokine concentration analysis in those subjects had been performed following the same methodology as in this study.

It should be noted that information about the concentration of some cytokines was missing in some individuals’ tear samples. Specifically, levels of IL-7 were missing in all the subjects.

The percentage of cytokine tear concentration values falling within the lower and upper RIs of each sample was determined.

3. Results

A total of 157 subjects were included in the study. The mean age of the sample was 46.36 ± 17.7 years (range, 18–88), and approximately half were women (n = 82; 52.2%; CI95%: 44.1%, 60.2%). There were no significant sex differences in age (men: 48.65 ± 19.19; women: 44.26 ± 16.04; Student’s t test p-value = 0.1202).

Individual tear samples were analyzed. The minimum detectable concentrations (in pg/mL) and detection rates for the molecules analyzed are shown in Table 1.

Eight (EGF, Fractalkine, IL-1RA, IL-7, IL-8, IP-10, MCP-1, and VEGF) out of the 18 molecules analyzed were detected in at least the 82% of the participants. Table 2 summarizes the distribution of the levels of these eight cytokines by age and sex.

The rest of the molecules analyzed (IL-1β, IL-2, IL-6, IL-10, IL-12p70, IL-15, IL-17A, IL-23, TNF-α, and IFN-γ) were not further considered for subsequent statistical analysis because their detection rates were below 50% (Table 1).

Table 3 shows the best GAMLSS models of each cytokine. Model selection involved the selection of the distribution for the corresponding cytokine and the selection of relevant explicative factors by the comparison of fitted models according to the smallest AIC. Four different distribution families were considered: Normal, Box-Cox Cole and Green, Power Exponential, and Box-Cox power exponential. For IL-8, IP-10, and VEGF, Power Exponential was chosen as the most appropriate distribution. Normal distribution was selected for Fractalkine, IL-7, and MCP-1. Finally, Box-Cox power exponential was the less selected distribution, chosen in two of the eight cytokines considered (EGF and IL-1RA).

The potential explanatory variables were age and sex. The age influence on the location parameter was common to EGF, IL-1RA, IL-7, IL-8, IP-10, MCP-1, and VEGF, but not for Fractalkine, whose mean parameter depends exclusively on sex. In addition, age was added to model variability of IP-10 and skewness of VEGF. Sex had a relevant influence on location parameter of Fractalkine, IL-7, IL-8, and IP-10 and on variability of IP-10 (Table 3). A residual analysis of each of the final fitted models was carried out. In this diagnostic stage, none of the cytokines showed inadequacies in the final fitted model (see Residual Analysis in Appendix A).

The age- and sex-related 90% reference (5%, 95%) and median (50%) curves calculated using the best-fitted GAMLSS model for each cytokine are shown in Figure 2; complete numerical values can be found in Appendix B. The levels of IL-1RA, IL-8, IP-10, and MCP-1 showed an age-related increase; in the case of IP-10 and MCP-1 this increase was approximately linear but for IL-1RA and IL-8, the increase was only relevant in older subjects (70 years old or more). In contrast, EGF tear concentration decreased linearly with age.

Age patterns were also observed for IL-7 and VEGF tear levels. IL-7 concentration increased up to approximately 40 years, then a mild decrease was observed up to 70 years and then increased again. For VEGF, higher levels were observed in the middle age range (30 to 60 years approximately).

Regarding sex influence, Fractalkine levels were generally higher in men, whereas IL-7, IL-8, and IP-10 tear concentrations were higher in women (Figure 2).

In order to validate the estimated RIs, we used them to evaluate the tear cytokine concentrations from an external independent sample from our laboratory database of 40 healthy individuals (validation sample). In addition, tear cytokine concentration values from 13 cases of severe DED (utility sample; also from our own database) were also used to show the utility of the estimated RIs for diagnosis of DED. RIs were evaluated for all cytokines except for IL-7 whose levels were missing in all subjects of both samples. Table 4 summarizes the proportion of values that fell within the appropriate age- and sex-specific RIs for all other cytokines in each sample. The complete demographic information and classification for each subject in the validation and utility samples is shown in Table A11 and Table A12 (respectively) in Appendix C. More than 92% of validation sample values fell within the appropriate age- and sex-specific RIs for each cytokine. Concerning the utility sample, IL-1RA and IL-8 levels were above the upper RIs in 92.3% of those patients diagnosed with severe DED compared to 7.7% and 5% of those clinically classified as not having DED. Fractalkine and VEGF tear levels were also increased in the DED sample; in these cases, the percentages of values over the upper RIs were 72.7% and 75% versus 5.6% and 6.3% of the validation sample, respectively. On the other hand, EGF and IP-10 showed lower concentrations in the utility sample, finding that 23.1% and 38.5% of the DED patients had levels below the corresponding lower RI. These percentages are significantly higher than those found in the validation sample (0% and 2.9%, respectively). Interestingly, all the patients of severe DED sample had tear cytokine concentrations that fell outside the RIs in at least two of the cytokines evaluated (see Table A12 in Appendix C).

4. Discussion

In this work, we have established the sex- and age-adjusted RIs of eight (EGF, Fractalkine, IL-1RA, IL-7, IL-8, IP-10, MCP-1 and VEGF) tear cytokine levels. RIs are used to describe the normal range of a parameter, and their important role in clinical medicine is beyond discussion as an initial indication for further investigation or treatment. They provide a basis for comparison of measured values of patients to a reference, and they could be very useful to diagnose the condition of a subject in epidemiological and clinical studies.

The most commonly used definition of a RI is the interval of values containing the central 90% or 95% of a healthy population [28]. It is possible to estimate population percentiles from an adequate sample by computing observed percentiles. As a result, single cut-off values will be obtained. However, this result may not be very useful in the case of cytokine levels, since the distribution of cytokine concentrations, as inflammatory marker measurements, could be heavily dependent on individual characteristics, such as age and sex among others [29]. In addition, immunological data very frequently have asymmetrical distributions [30] and statistical techniques that take into account this issue should be chosen. For these reasons, in this study, GAMLSS models [17] were applied to derive age- and sex-dependent centile curves. The main advantage of these type of models is that the shape of the distribution, rather than just the mean or the variance, for the response variable can vary according to more than one explanatory variable. In our case, for example, the skewness parameter for VEGF concentration varied according to the subject age and it would be wrongly modeled by other more classical statistical methodology.

Our results showed that tear levels of all studied cytokines, except Fractalkine, presented an age-dependent behavior in any of the distributional parameters. Age-associated changes in human tear film composition and proteome have already been shown in other studies [12,15,16] in agreement with that reported as well in other fluids, such as plasma or cerebrospinal fluid [31]. Indeed, the relationship between the aging process and immune response is not new. Even, the expression “inflamm-aging” has been coined to describe the progressive increase in pro-inflammatory status as age increases [32,33]. There is a lot of literature about age, the role of inflammation markers and its link with the development of disease (reviewed in [31]). The pro- and anti-inflammatory cytokine balance changes with increasing age, and low-grade inflammation is more common among elderly population, leading to a chronic low-grade inflammatory phenotype [34]. Accordingly, age-dependent behavior of several cytokines levels has widely been reported in different biological fluids, such as plasma, cerebrospinal fluid, and also in tears, although tear studies are still scarce. Nevertheless, cytokine levels in these different biological sources have shown to have similar evolution patterns and common to our results in tears. For instance, in our study MCP-1 and IP-10 tear levels showed a sustained age-dependent increase, similar to what has already been observed in serum levels [35]. In addition, IP-10 increased linearly with age in cerebrospinal fluid [36], and plasma MCP-1 levels were significantly higher in an elderly group (70 to 92 years) compared to adults (21 to 50 years) [37]. This latest study also showed that plasma levels of EGF were significantly lower in the elderly group, which is compatible with our results in tears. In our study, IL-1RA and IL-8 tear levels showed as well a very important increase but only in old ages (68 to 90 years). This same pattern has already been reported in tear [12] and in serum levels [38] of IL-8, and in plasma levels of IL-1RA [39]. Even, these two cytokines have been proposed as markers linked to longevity in very old people [38,40]. An increase of IL-8 levels with age has recently been described in cerebrospinal fluid as well, but, in this case, the evolution is practically linear [36].

The role of IL-7 in immunosenescence, the aging of the immune system, has been known for long [41]. In this work, we have found a cyclic evolution of IL-7 tear levels with higher levels around 30 years old and increasing again after 70 years old. This result disagrees with the findings by Micera et al. [12] which found an increase of IL-7 in middle age people (41–60 years old) with respect to young (18–40) and old (<60 years old). However, the observed behavior for IL-7 tear levels in our study could be consistent with the findings of Nasi et al. [42], that evaluated plasma IL-7 levels in three age groups: young (aged 20–45), middle-aged (aged 58–62), and centenarians (range 98–100 years). Despite the fact that they found no significant differences between the age groups, the evolution pattern of IL-7 levels was similar to ours. There is more confusion regarding VEGF evolution. In serum, for example, both an increase [43] and a decrease [44] of VEGF levels with aging have been found. In plasma samples of patients from 43 to 80 years old with stable coronary artery disease a significant negative correlation was observed [45]; in contrast, recently, Chakraborty et al. [46] detected a positive correlation in cerebrospinal fluid from subjects aged from 40 years. Other authors such as Rübenhagen et al. [47] found no association between the concentration of VEGF in synovial fluid and age. In our case, the relation between VEGF tear level and age can be described as an inverted-U shape curve, with maximum levels between 30–60 years old.

Our results regarding the age-dependence in tear cytokine levels are also in agreement with the findings by Nättinen et al. [16] for tear proteome in which the majority of the identified proteins in their study had the most notable increase/decrease in subjects aged 60 or more.

Having age-dependent RI of tear cytokine levels can be very useful for the diagnosis of ocular disorders. Several very common eye problems are related to the phenomenon of inflammaging. For example, an association between age and DED has been assumed for a long time, and recently, Di Zazzo et al. [48] have directly related inflammaging to a progressive impairment of ocular surface system function. An ocular disease in which age is an important risk factor is age-related macular degeneration, affecting people aged 65 years and older. Several inflammatory gene polymorphisms have been linked to this disorder indicating an important role of inflammation and immune-mediated processes, and, therefore, directly related to inflammaging term [49]. Finally, low-grade inflammation has been described as a causal factor in the pathogenesis of glaucoma [50]. Although glaucoma can affect people of all ages, is most common in adults in their 70s and 80s, so it is accepted that glaucoma risk increases with age [51].

Sex also has an influence on tear levels of some of the analyzed cytokines. Particularly, Fractalkine tear levels were higher in males compared to females and IL-7, IL-8, and IP-10 higher in females compared to males; whereas, EGF, IL-1RA, MCP-1, and VEGF levels were not affected by sex. Some other studies have also shown sex-dependence in cytokine levels in other fluids. For instance, Nasi et al. (2006) [42], found a sex-dependence for IL-7 serum levels, with higher levels in centenarian females than in males, in agreement with our results in tears. Moreover, the study by Nättinen et al. (2019) [16] revealed that protein-age correlations in tear film were, in most cases, more significant in males, although mammaglobin-A (SCGB2A2) protein had a higher correlation coefficient for females. Changes in sex hormones have been suggested as connected with these differences. In addition, very recently, stronger relationships between age and IL-6 and high-sensitivity C-reactive protein plasma levels have been found in females compared to males [52]. Although the cytokines analyzed in that work are not the same that those in this study, our results reflected the same behavior in relation to IL-7, IL-8, and IP-10 pro-inflammatory cytokines.

The association between tear concentration of cytokines and DED has been widely reported (reviewed in [3,4,53,54]). However, it is often difficult to use cytokines as diagnostic tools due to the lack of diagnostic cut-offs. In this study, we have shown the utility of our RIs to diagnose DED. All of the severe DED patients in our sample showed abnormal tear levels of at least two of the evaluated cytokines when compared to the established RIs. Fractalkine, IL-1RA, IL-8, and VEGF showed high levels in more than 70% of the DED patients. Moreover, EGF and IP-10 were found to be infra-expressed in 23.1% and 38.5% of patients respectively of our DED sample. Altered tear levels of all these molecules have already been reported in different clinical subgroups of patients with DED. So, for example, tear concentrations of Fractalkine, IL-1RA, IL-8, and VEGF have frequently been found significantly increased in DED patients as compared to normal controls [8,18,55,56,57,58,59], whereas tear EGF level has been found significantly decreased in the more severe forms of DED [55,60,61]. IP-10 tear levels were also found decreased compared to normal in 38.5% of the DED patients’ sample, in agreement with results in tears from other studies involving severe cases of DED associated to ocular graft versus host disease [62], Stevens-Johnson syndrome [63], and, as preliminary findings suggest, in primary Sjögren’s syndrome subjects [64,65]. However, several authors have also found tear IP-10 concentration to be overexpressed in tears of patients with moderate evaporative DED [21] and with Sjögren’s syndrome [66]. It has been suggested that the anti-angiogenic properties of tear IP-10 may contribute to normal ocular surface immune privilege [65] and the low levels in Stevens-Johnson syndrome have been suggested to lead to a severe smoldering inflammation and progressive subconjunctival fibrosis in these patients [63]. Further studies with a bigger DED sample are warranted to evaluate the role of IP-10 in DED.

A critical issue of establishing RIs is to recruit a valid group of reference individuals. Two fundamental aspects are important in this regard: the definition of “reference individual” and the sample size. In this case, a reference individual should be a subject without eye pathologies. However, the characterization of this type of subject is not always easy, so collecting medium/large samples from the reference population can be quite expensive. Certainly, although above 150 individuals, our sample could be considered not sufficiently large to properly represent the reference population. Nevertheless, it is very important to highlight that a single ophthalmologist evaluated all subjects included in our sample, so the eye health criterion was homogeneous. Moreover, sample-handling factors such as tear collection, storage, or processing protocol may be a strong influence on the value of cytokine level and they could help to invalidate the sample as a reference sample. In this study, the tear samples were collected and processed by a single person, so this source of variability was also controlled. It also should be considered that in this study environmental conditions during subjects’ clinical evaluation and sample collection were not controlled; it is well known the fact that ocular surface parameters and tear mediator concentration are altered under adverse environmental conditions, particularly low humidity and/or high temperatures. Nevertheless, all subjects were evaluated and samples were obtained in the same air-conditioned rooms (in a normal range of temperature and humidity), assuring that all of them were evaluated under the same conditions. In general, aware that it will not always be possible to collect samples under controlled environmental conditions, our recommendation is to at least do it in the same office, under similar ambient conditions.

Additionally, it is noted that our validation, based on an external sample of subjects without ocular pathology, was appropriate for all of studied cytokines. As per clinical and laboratory standards institute C28-A3 guidelines [28] validation can be considered successful, since more than 20 samples were evaluated and more than 90% of values fell within the RIs. In some way, this result reinforces the validity of our reference sample. Finally, although this study is related to the use of specific cytokine/chemokine X-MAP-Luminex kits for molecule determination, the results obtained about the dependence on sex and age in tear cytokine levels would remain the same regardless of the use of kits from other vendors or another molecule determination assay.

5. Conclusions

We have estimated adjusted RIs of eight tear cytokine levels by GAMLSS regression models. We were able to determine the effect of age and sex in the concentration of each of these tear molecules in subjects without ocular disorders and adjusted their cut-off values taking into account these dependencies. Finally, we have showed the effectiveness of our RIs as diagnostic tools for severe DED. These age- and sex-intervals for tear cytokine levels could therefore be used as reference values in other tear cytokine studies in patients, and can be an important tool to extend the use of tear cytokines as biomarkers of ocular diseases.

Author Contributions

Study design: I.F., A.E.-d.-S., M.C., J.M.H.; Conduct of the study: A.P., A.E.-d.-S., C.G.-V., I.F.; Analysis: I.F., A.E.-d.-S.; Writing—original draft preparation: I.F., A.E.-d.-S., A.P.; Writing—review and editing: I.F., A.E.-d.-S., M.C., J.M.H.; Final manuscript approval: I.F., A.E.-d.-S., A.P., C.G.-V., M.C., J.M.H. All authors have read and agreed to the published version of the manuscript.

Funding

This research received no external funding.

Institutional Review Board Statement

The study was approved by the Institute of Applied Ophthalmobiology (IOBA) institutional review board and by the Ethics Committee of the University of Valladolid (16 January 2009) and complied with the Tenets of the Declaration of Helsinki.

Informed Consent Statement

Informed consent was obtained from all subjects involved in the study.

Data Availability Statement

Data is contained within the article.

Conflicts of Interest

The authors declare no conflict of interest. The funders had no role in the design of the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, or in the decision to publish the results.

Appendix A. Residual Analysis

This appendix shows the diagnostic tools based on normalized quantile residuals [25] used for checking the adequacy of each one of the fitted generalized additive models for location, scale and shape (GAMLSS) models. For a particular model, the following diagnostic tools are used:

Depending on the number of explanatory variables, four or five residuals plots for checking the normalized quantile residuals of the fitted model:
-
Histogram of the residuals adding their kernel smooth density estimate (dotted line) and standard normal density (solid line);
-
QQ-normal plot of the residuals;
-
Residuals against the fitted values of the μ parameter;
-
Residuals against each one of the included explanatory variables.

If the residuals behave well, the density estimate of the residuals will be approximately normal, the normal Q-Q plot will be approximately linear (with intercept 0 and gradient 1) and the bottom two or three plots will show a random scatter around the horizontal line at 0, and similar distribution of residuals for both sexes. In these plots line on Y = 0 (solid line) and fitted smooth curve (dotted line) are also represented.

Summary statistics of the normalized quantile residuals are also calculated to check deviation from normally assumption. Specifically, mean, variance, skewness, and kurtosis coefficients are estimated and theirs reference values are 0, 1, 0, and 3, respectively.
Worm plots [26] of the normalized quantile residuals to identify adjustments for some of fitted model parameters. Below are detailed the main features:
-
Points on the plot show the gap between the ordered residuals and their approximate expected values from a standard normal distribution, represented by the horizontal dotted line.
-
Two elliptic curves represent the approximate point-wise 95% confidence bands. When the model is correct, approximately 95% of the points to lie between the confidence bands. A higher percentage of the points outside the elliptic curves indicates that the fitted terms are inadequate.
-
The shape of cubic fitted curve to the points (solid curve) shows inadequacies in the model described in Table A1.

Multiple worm plots are also built in order to identify failures of the model within different ranges of the explanatory variables. Age is cut into four non-overlapping intervals with similar number of observations and sex splits the data into female and male groups.

Table A1. Deficiencies in the residuals and in the fitted distribution according to the shapes for the worm plot of the normalized quantile residuals.

Shape of the Fitted Curve in the Worm Plot	Deficiencies in the Residuals	Deficiencies in the Fitted Distribution
Level above the origin	Mean too high	Fitted location parameter (μ) is too low
Level below the origin	Mean too low	Fitted location parameter (μ) is too high
Line with positive slope	Variance too high	Fitted scale parameter (σ) is too low
Line with negative slope	Variance too low	Fitted scale parameter (σ) is too high
Curve with U shape	Positive skewness	Fitted skewness parameter (ν) is too low
Curve with inverted U shape	Negative skewness	Fitted skewness parameter (ν) is too high
Curve with S shape with left bent down	Leptokurtosis (more peaked than that of a normal distribution)	Fitted kurtosis parameter (τ) is too low
Curve with S shape with left bent up	Platykurtosis (flatter than that of a normal distribution)	Fitted kurtosis parameter (τ) is too high

Table A2. Summary statistics of the quantile residuals for each fitted model. In all cases, the statistics support that the residuals are approximately normally distributed: the mean is nearly zero, the variance nearly one, the coefficient of skewness is near zero and the coefficient of kurtosis is near three.

	Mean	Variance	Skewness Coef.	Kurtosis Coef.
	(95% CI)	(95% CI)	(95% CI)	(95% CI)
EGF	0	1	−0.01	2.85
EGF	(−0.156, 0.16)	(0.814, 1.27)	(−0.324, 0.286)	(2.398, 3.315)
Fractalkine	0	1.01	0.05	2.74
Fractalkine	(−0.158, 0.158)	(0.816, 1.273)	(−0.269, 0.33)	(2.392, 3.15)
IL-1RA	0	1.01	0	2.92
IL-1RA	(−0.158, 0.159)	(0.816, 1.274)	(−0.3, 0.325)	(2.416, 3.492)
IL-7	0	1.01	−0.2	2.75
IL-7	(−0.158, 0.158)	(0.816, 1.273)	(−0.543, 0.117)	(2.284, 3.207)
IL-8	0.04	1	0.16	2.89
IL-8	(−0.113, 0.203)	(0.814, 1.271)	(−0.219, 0.453)	(2.409, 3.399)
IP-10	0.05	1	−0.15	2.97
IP-10	(−0.108, 0.208)	(0.813, 1.269)	(−0.442, 0.215)	(2.448, 3.496)
MCP-1	0	1.01	0.08	2.68
MCP-1	(−0.158, 0.158)	(0.816, 1.273)	(−0.215, 0.375)	(2.269, 3.16)
VEGF	0.04	1.05	−0.07	2.83
VEGF	(−0.121, 0.203)	(0.855, 1.334)	(−0.369, 0.229)	(2.42, 3.288)

CI = Confidence interval; EGF = Epidermal Growth Factor; IL = Interleukin; RA = Receptor Antagonist; IP = Induced Protein; MCP = Monocyte Chemotactic Protein; VEGF = Vascular Endothelial Growth Factor.

Figure A1. Residual plots from the fitted model to EGF levels. The normalized quantile residuals behave well. The histogram and normal QQ-plot (two top plots) do not show deviates from normality assumption. Residuals against fitted values and against age plots (two bottom plots) show a random scatter around the horizontal line at 0.

Figure A2. Worm plots from the fitted model to EGF levels. General worm plot is on the top of the figure. Age is included in the model, and four individual worm plots, corresponding to the four age intervals, are shown along rows from top left to bottom right. The five worm plots suggest an adequate model because the points in each of the plots lie within the two elliptic confidence bands.

Figure A3. Residual plots from the fitted model to Fractalkine levels. The normalized quantile residuals behave well. The histogram and normal QQ-plot (two top plots) do not show deviates from normality assumption. Residuals against fitted values (left bottom plot) shows a random scatter around the horizontal line at 0. Boxplot that compares the residuals for the two sex levels (right bottom plot) shows similar distributions in both male and female.

Figure A4. Worm plots from the fitted model to Fractalkine levels. General worm plot is on the top of the figure. Sex is included in the model, and two individual worm plots, corresponding to the female and male groups, are shown on the bottom on the figure. The three worm plots suggest an adequate model because the points in each of the plots lie within the two elliptic confidence bands.

Figure A5. Residual plots from the fitted model to IL-1RA levels. The normalized quantile residuals behave well. The histogram and normal QQ-plot (two top plots) do not show deviates from normality assumption. Residuals against fitted values (left bottom plot) and against age (center bottom plot) show a random scatter around the horizontal line at 0.

Figure A6. Worm plots from the fitted model to IL-1RA levels. General worm plot is on the top of the figure. Age is included in the model, and four individual worm plots, corresponding to the four age intervals, are shown on the bottom on the figure, respectively. The five worm plots suggest an adequate model because the points in each of the plots lie within the two elliptic confidence bands.

Figure A7. Residual plots from the fitted model to IL-7 levels. The normalized quantile residuals behave well. The histogram and normal QQ-plot (two top plots) do not show deviates from normality assumption. Residuals against fitted values (left bottom plot) and against age (center bottom plot) show a random scatter around the horizontal line at 0. Boxplot that compares the residuals for the two sex levels (right bottom plot) shows similar distributions in both male and female.

Figure A8. Worm plots from the fitted model to IL-7 levels. General worm plot is on the top of the figure. Age and sex are included in the model, and four and two individual worm plots, corresponding to the four age intervals (left) and to the female and male groups (right), are shown on the bottom on the figure, respectively. The seven worm plots suggest an adequate model because the points in each of the plots lie within the two elliptic confidence bands.

Figure A9. Residual plots from the fitted model to IL-8 levels. The normalized quantile residuals behave well. The histogram and normal QQ-plot (two top plots) do not show deviates from normality assumption. Residuals against fitted values (left bottom plot) and against age (center bottom plot) show a random scatter around the horizontal line at 0. Boxplot that compares the residuals for the two sex levels (right bottom plot) shows similar distributions in both male and female.

Figure A10. Worm plots from the fitted model to IL-8 levels. General worm plot is on the top of the figure. Age and sex are included in the model, and four and two individual worm plots, corresponding to the four age intervals (left) and to the female and male groups (right), are shown on the bottom on the figure, respectively. The seven worm plots suggest an adequate model because the points in each of the plots lie within the two elliptic confidence bands.

Figure A11. Residual plots from the fitted model to IP-10 levels. The normalized quantile residuals behave well. The histogram and normal QQ-plot (two top plots) do not show deviates from normality assumption. Residuals against fitted values (left bottom plot) and against age (center bottom plot) show a random scatter around the horizontal line at 0. Boxplot that compares the residuals for the two sex levels (right bottom plot) shows similar distributions in both male and female.

Figure A12. Worm plots from the fitted model to IP-10 levels. General worm plot is on the top of the figure. Age and sex are included in the model, and four and two individual worm plots, corresponding to the four age intervals (left) and to the female and male groups (right), are shown on the bottom on the figure, respectively. The seven worm plots suggest an adequate model because the points in each of the plots lie within the two elliptic confidence bands.

Figure A13. Residual plots from the fitted model to MCP-1 levels. The normalized quantile residuals behave well. The histogram and normal QQ-plot (two top plots) do not show deviates from normality assumption. Residuals against fitted values and against age plots (two bottom plots) show a random scatter around the horizontal line at 0.

Figure A14. Worm plots from the fitted model to MCP-1 levels. General worm plot is on the top of the figure. Age is included in the model, and four individual worm plots, corresponding to the four age intervals, are shown along rows from top left to bottom right. The five worm plots suggest an adequate model because the points in each of the plots lie within the two elliptic confidence bands.

Figure A15. Residual plots from the fitted model to VEGF levels. The normalized quantile residuals behave well. The histogram and normal QQ-plot (two top plots) do not show deviates from normality assumption. Residuals against fitted values and against age plots (two bottom plots) show a random scatter around the horizontal line at 0.

Figure A16. Worm plots from the fitted model to VEGF levels. General worm plot is on the top of the figure. Age is included in the model, and four individual worm plots, corresponding to the four age intervals, are shown along rows from top left to bottom right. The five worm plots suggest an adequate model because the points in each of the plots lie within the two elliptic confidence bands.

Appendix B. Centile Curves Estimation

This appendix shows the estimations of centile curves for each evaluated cytokine based on the fitted generalized additive models for location, scale and shape (GAMLSS) models.

Table A3. Percentiles of EGF levels estimated with the fitted GAMLSS model. Family distribution Box-Cox Power Exponential with age linear dependence for location parameter (μ).

Age (Years)	Percentile
Age (Years)	2.5th	5th	50th	95th	97.5th
18	32.55	50.26	508.25	3314.48	4333.62
19	32.32	49.86	501.91	3260.74	4261.06
20	32.10	49.47	495.64	3207.88	4189.71
21	31.87	49.08	489.45	3155.87	4119.56
22	31.65	48.70	483.34	3104.71	4050.58
23	31.43	48.31	477.31	3054.37	3982.76
24	31.21	47.93	471.35	3004.86	3916.07
25	30.99	47.56	465.46	2956.14	3850.50
26	30.77	47.18	459.65	2908.21	3786.03
27	30.56	46.81	453.91	2861.07	3722.63
28	30.34	46.44	448.24	2814.68	3660.30
29	30.13	46.08	442.65	2769.05	3599.01
30	29.92	45.71	437.12	2724.16	3538.75
31	29.71	45.35	431.66	2679.99	3479.50
32	29.50	45.00	426.27	2636.54	3421.24
33	29.30	44.64	420.95	2593.80	3363.95
34	29.09	44.29	415.70	2551.75	3307.63
35	28.89	43.94	410.51	2510.38	3252.24
36	28.69	43.60	405.38	2469.68	3197.79
37	28.49	43.26	400.32	2429.64	3144.25
38	28.29	42.92	395.32	2390.25	3091.60
39	28.09	42.58	390.39	2351.50	3039.83
40	27.89	42.24	385.51	2313.38	2988.93
41	27.70	41.91	380.70	2275.87	2938.89
42	27.50	41.58	375.94	2238.97	2889.68
43	27.31	41.25	371.25	2202.67	2841.29
44	27.12	40.93	366.62	2166.96	2793.72
45	26.93	40.61	362.04	2131.83	2746.94
46	26.74	40.29	357.52	2097.27	2700.95
47	26.55	39.97	353.05	2063.27	2655.72
48	26.37	39.66	348.65	2029.82	2611.26
49	26.18	39.34	344.29	1996.91	2567.53
50	26.00	39.03	339.99	1964.54	2524.54
51	25.82	38.73	335.75	1932.69	2482.27
52	25.64	38.42	331.56	1901.35	2440.71
53	25.46	38.12	327.42	1870.53	2399.84
54	25.28	37.82	323.33	1840.20	2359.66
55	25.10	37.52	319.29	1810.37	2320.15
56	24.93	37.23	315.31	1781.02	2281.30
57	24.75	36.93	311.37	1752.14	2243.10
58	24.58	36.64	307.48	1723.74	2205.54
59	24.41	36.36	303.64	1695.79	2168.61
60	24.24	36.07	299.85	1668.30	2132.30
61	24.07	35.79	296.11	1641.25	2096.60
62	23.90	35.50	292.41	1614.65	2061.49
63	23.73	35.22	288.76	1588.47	2026.98
64	23.57	34.95	285.15	1562.72	1993.04
65	23.40	34.67	281.59	1537.38	1959.67
66	23.24	34.40	278.08	1512.46	1926.85
67	23.07	34.13	274.61	1487.94	1894.59
68	22.91	33.86	271.18	1463.81	1862.87
69	22.75	33.59	267.79	1440.08	1831.68
70	22.59	33.33	264.45	1416.73	1801.01
71	22.44	33.07	261.15	1393.77	1770.85
72	22.28	32.81	257.89	1371.17	1741.20
73	22.12	32.55	254.67	1348.94	1712.04
74	21.97	32.29	251.49	1327.07	1683.38
75	21.81	32.04	248.35	1305.56	1655.19
76	21.66	31.79	245.25	1284.39	1627.48
77	21.51	31.54	242.18	1263.57	1600.23
78	21.36	31.29	239.16	1243.08	1573.43
79	21.21	31.04	236.17	1222.93	1547.09
80	21.06	30.80	233.22	1203.10	1521.18
81	20.91	30.56	230.31	1183.60	1495.71
82	20.77	30.32	227.44	1164.41	1470.67
83	20.62	30.08	224.60	1145.53	1446.04
84	20.48	29.84	221.79	1126.96	1421.83
85	20.33	29.61	219.02	1108.69	1398.02
86	20.19	29.37	216.29	1090.71	1374.61
87	20.05	29.14	213.59	1073.03	1351.60
88	19.91	28.91	210.92	1055.63	1328.97

Table A4. Percentiles of Fractalkine levels estimated with the fitted GAMLSS model. Family distribution Normal with sex dependence for location parameter (μ).

Percentile for Women					Percentile for Men
2.5th	5th	50th	95th	97.5th	2.5th	5th	50th	95th	97.5th
19.25	26.64	145.21	791.44	1095.20	25.93	35.88	195.53	1065.69	1474.72

Table A5. Percentiles of IL-1RA levels estimated with the fitted GAMLSS model. Family distribution Box-Cox Power Exponential with age dependence modelled as a cubic spline for location parameter (μ).

Age (Years)	Percentile
Age (Years)	2.5th	5th	50th	95th	97.5th
18	21.92	29.06	238.21	4675.57	8526.01
19	22.51	29.92	249.70	5028.39	9216.92
20	23.11	30.78	261.59	5402.87	9954.05
21	23.71	31.66	273.72	5794.25	10,728.33
22	24.30	32.52	285.90	6197.32	11,529.63
23	24.88	33.36	298.08	6609.48	12,352.84
24	25.42	34.16	309.71	7011.97	13,160.26
25	25.89	34.84	319.89	7370.82	13,882.94
26	26.24	35.35	327.47	7642.29	14,431.31
27	26.42	35.62	331.57	7790.47	14,731.23
28	26.46	35.67	332.35	7819.00	14,789.03
29	26.37	35.55	330.52	7752.64	14,654.63
30	26.21	35.31	326.90	7621.70	14,389.67
31	26.00	35.01	322.35	7458.50	14,059.90
32	25.79	34.69	317.64	7290.98	13,721.93
33	25.60	34.41	313.52	7145.39	13,428.64
34	25.46	34.22	310.57	7042.10	13,220.84
35	25.38	34.10	308.87	6982.50	13,101.03
36	25.34	34.04	307.94	6950.29	13,036.30
37	25.31	33.99	307.33	6928.79	12,993.11
38	25.29	33.97	306.89	6913.76	12,962.94
39	25.28	33.95	306.63	6904.68	12,944.69
40	25.28	33.94	306.57	6902.55	12,940.43
41	25.28	33.94	306.56	6902.30	12,939.93
42	25.27	33.93	306.39	6896.08	12,927.42
43	25.24	33.89	305.73	6873.36	12,881.82
44	25.18	33.80	304.45	6829.08	12,792.95
45	25.09	33.67	302.54	6762.95	12,660.31
46	24.97	33.49	299.97	6674.32	12,482.68
47	24.81	33.26	296.54	6557.01	12,247.82
48	24.61	32.97	292.43	6417.22	11,968.36
49	24.40	32.66	287.91	6264.80	11,664.14
50	24.17	32.33	283.26	6109.30	11,354.31
51	23.96	32.01	278.72	5958.78	11,054.94
52	23.75	31.70	274.40	5816.63	10,772.72
53	23.54	31.40	270.20	5679.81	10,501.53
54	23.33	31.09	265.87	5539.90	10,224.70
55	23.09	30.75	261.15	5388.59	9925.87
56	22.82	30.36	255.77	5218.18	9590.03
57	22.51	29.92	249.70	5028.20	9216.57
58	22.17	29.41	242.94	4819.77	8807.96
59	21.79	28.86	235.60	4596.83	8372.30
60	21.39	28.30	228.11	4373.04	7936.51
61	21.01	27.74	220.93	4162.42	7527.78
62	20.67	27.24	214.51	3977.06	7169.31
63	20.38	26.83	209.24	3827.22	6880.39
64	20.17	26.54	205.54	3723.27	6680.42
65	20.07	26.39	203.73	3672.81	6583.48
66	20.08	26.41	203.92	3678.21	6593.86
67	20.21	26.59	206.21	3742.14	6716.69
68	20.46	26.95	210.79	3871.14	6964.99
69	20.86	27.52	218.05	4078.99	7366.30
70	21.40	28.31	228.28	4378.15	7946.44
71	22.09	29.30	241.42	4773.25	8716.92
72	22.90	30.48	257.36	5268.45	9689.03
73	23.82	31.82	276.02	5869.64	10,877.90
74	24.85	33.32	297.48	6589.20	12,312.25
75	25.98	34.98	321.91	7442.84	14,028.28
76	27.22	36.79	349.47	8449.34	16,069.59
77	28.55	38.76	380.35	9629.55	18,485.24
78	29.98	40.88	414.72	11,005.47	21,328.07
79	31.49	43.15	452.69	12,599.45	24,653.07
80	33.10	45.55	494.35	14,433.69	28,516.30
81	34.78	48.09	539.83	16,534.26	32,983.42
82	36.55	50.76	589.47	18,939.52	38,148.14
83	38.41	53.59	643.71	21,696.71	44,126.04
84	40.37	56.58	703.05	24,861.31	51,053.90
85	42.44	59.75	768.06	28,498.56	59,094.05
86	44.62	63.10	839.37	32,685.48	68,439.71
87	46.92	66.66	917.71	37,513.16	79,321.44
88	49.36	70.45	1003.80	43,082.66	91,999.10

Table A6. Percentiles of IL-7 levels estimated with the fitted GAMLSS model. Family distribution Power Exponential with age dependence modelled as a cubic spline for location and skewness parameters (μ and ν respectively) and as linear for scale parameter (σ). Sex is added as a qualitative variable for location and skewness parameters (μ and ν respectively).

Age (Years)	Percentile for Women					Percentile for Men
Age (Years)	2.5th	5th	50th	95th	97.5th	2.5th	5th	50th	95th	97.5th
18	11.48	13.94	38.35	105.53	128.11	8.88	10.79	29.68	81.65	99.12
19	11.86	14.40	39.61	109.00	132.32	9.18	11.14	30.65	84.33	102.38
20	12.25	14.87	40.92	112.59	136.68	9.48	11.51	31.66	87.12	105.76
21	12.66	15.37	42.28	116.32	141.21	9.79	11.89	32.71	90.00	109.26
22	13.08	15.88	43.68	120.19	145.91	10.12	12.28	33.80	93.00	112.90
23	13.52	16.41	45.15	124.22	150.80	10.46	12.70	34.93	96.11	116.68
24	13.96	16.95	46.64	128.33	155.79	10.80	13.12	36.09	99.30	120.55
25	14.41	17.50	48.14	132.45	160.79	11.15	13.54	37.25	102.48	124.41
26	14.85	18.02	49.59	136.45	165.65	11.49	13.95	38.37	105.58	128.17
27	15.26	18.52	50.96	140.21	170.21	11.80	14.33	39.43	108.49	131.70
28	15.63	18.97	52.20	143.63	174.36	12.09	14.68	40.39	111.13	134.91
29	15.95	19.36	53.28	146.59	177.96	12.34	14.98	41.22	113.43	137.69
30	16.22	19.69	54.17	149.04	180.93	12.55	15.23	41.91	115.32	139.99
31	16.42	19.94	54.86	150.94	183.23	12.71	15.43	42.45	116.79	141.77
32	16.57	20.11	55.34	152.27	184.85	12.82	15.56	42.82	117.82	143.03
33	16.65	20.22	55.62	153.04	185.78	12.88	15.64	43.04	118.41	143.75
34	16.67	20.24	55.69	153.24	186.03	12.90	15.66	43.09	118.57	143.94
35	16.64	20.20	55.58	152.94	185.66	12.88	15.63	43.01	118.33	143.65
36	16.56	20.11	55.32	152.22	184.78	12.82	15.56	42.81	117.78	142.98
37	16.45	19.97	54.93	151.15	183.49	12.73	15.45	42.51	116.95	141.97
38	16.30	19.79	54.45	149.80	181.86	12.61	15.31	42.13	115.91	140.71
39	16.13	19.58	53.88	148.24	179.95	12.48	15.15	41.69	114.70	139.24
40	15.94	19.35	53.24	146.50	177.84	12.33	14.97	41.20	113.35	137.61
41	15.73	19.10	52.55	144.60	175.54	12.17	14.78	40.66	111.88	135.82
42	15.51	18.83	51.80	142.53	173.03	12.00	14.57	40.08	110.28	133.88
43	15.26	18.53	50.98	140.27	170.29	11.81	14.34	39.45	108.54	131.76
44	15.00	18.21	50.10	137.84	167.33	11.60	14.09	38.76	106.65	129.47
45	14.72	17.87	49.17	135.30	164.25	11.39	13.83	38.05	104.68	127.08
46	14.44	17.53	48.23	132.70	161.10	11.17	13.56	37.32	102.68	124.65
47	14.16	17.18	47.28	130.10	157.93	10.95	13.30	36.59	100.66	122.20
48	13.88	16.85	46.35	127.53	154.81	10.74	13.03	35.86	98.67	119.79
49	13.60	16.51	45.44	125.02	151.77	10.53	12.78	35.16	96.73	117.43
50	13.34	16.19	44.54	122.55	148.78	10.32	12.53	34.46	94.83	115.11
51	13.07	15.87	43.66	120.14	145.84	10.11	12.28	33.78	92.95	112.84
52	12.82	15.56	42.82	117.82	143.03	9.92	12.04	33.13	91.16	110.67
53	12.58	15.28	42.04	115.66	140.40	9.74	11.82	32.52	89.49	108.64
54	12.37	15.01	41.31	113.65	137.97	9.57	11.62	31.96	87.94	106.75
55	12.17	14.77	40.64	111.81	135.74	9.41	11.43	31.44	86.51	105.03
56	11.98	14.55	40.02	110.12	133.68	9.27	11.25	30.97	85.20	103.43
57	11.81	14.34	39.44	108.52	131.74	9.14	11.09	30.52	83.97	101.93
58	11.64	14.14	38.89	107.01	129.91	9.01	10.94	30.09	82.80	100.52
59	11.49	13.95	38.37	105.58	128.17	8.89	10.79	29.69	81.69	99.17
60	11.34	13.77	37.88	104.22	126.52	8.77	10.65	29.31	80.64	97.89
61	11.20	13.60	37.41	102.93	124.96	8.67	10.52	28.95	79.64	96.69
62	11.07	13.44	36.97	101.73	123.50	8.56	10.40	28.61	78.71	95.55
63	10.95	13.29	36.57	100.62	122.15	8.47	10.28	28.30	77.85	94.51
64	10.84	13.16	36.21	99.62	120.93	8.39	10.18	28.01	77.08	93.57
65	10.74	13.04	35.89	98.74	119.86	8.31	10.09	27.77	76.40	92.74
66	10.66	12.94	35.61	97.98	118.95	8.25	10.01	27.55	75.81	92.03
67	10.59	12.86	35.38	97.34	118.17	8.20	9.95	27.37	75.32	91.43
68	10.54	12.79	35.19	96.83	117.55	8.15	9.90	27.23	74.92	90.95
69	10.50	12.74	35.06	96.47	117.11	8.12	9.86	27.13	74.64	90.61
70	10.47	12.72	34.99	96.27	116.86	8.10	9.84	27.07	74.49	90.42
71	10.47	12.71	34.98	96.24	116.83	8.10	9.84	27.06	74.46	90.40
72	10.49	12.73	35.03	96.39	117.02	8.12	9.85	27.11	74.58	90.54
73	10.53	12.78	35.16	96.73	117.43	8.14	9.89	27.20	74.85	90.86
74	10.59	12.85	35.36	97.29	118.11	8.19	9.94	27.36	75.28	91.38
75	10.67	12.96	35.65	98.10	119.08	8.26	10.03	27.59	75.90	92.14
76	10.79	13.10	36.04	99.16	120.37	8.35	10.13	27.88	76.72	93.14
77	10.93	13.27	36.52	100.48	121.98	8.46	10.27	28.26	77.74	94.38
78	11.11	13.48	37.09	102.06	123.90	8.59	10.43	28.70	78.97	95.87
79	11.31	13.73	37.77	103.92	126.16	8.75	10.62	29.22	80.41	97.61
80	11.54	14.01	38.54	106.05	128.74	8.93	10.84	29.82	82.06	99.62
81	11.80	14.33	39.42	108.45	131.66	9.13	11.08	30.50	83.91	101.87
82	12.09	14.68	40.38	111.10	134.87	9.35	11.36	31.24	85.96	104.36
83	12.40	15.06	41.43	113.99	138.38	9.60	11.65	32.05	88.20	107.07
84	12.74	15.47	42.56	117.09	142.15	9.86	11.97	32.93	90.60	109.98
85	13.10	15.90	43.76	120.40	146.16	10.14	12.31	33.86	93.16	113.09
86	13.48	16.36	45.02	123.88	150.39	10.43	12.66	34.84	95.85	116.36
87	13.88	16.84	46.35	127.52	154.80	10.74	13.03	35.86	98.66	119.78
88	14.28	17.34	47.71	131.28	159.37	11.05	13.42	36.92	101.58	123.31

Table A7. Percentiles of IL-8 levels estimated with the fitted GAMLSS model. Family distribution Power Exponential with age dependence modelled as a cubic spline and sex as a qualitative variable for location parameter (μ).

Age (Years)	Percentile for Women					Percentile for Men
Age (Years)	2.5th	5th	50th	95th	97.5th	2.5th	5th	50th	95th	97.5th
18	16.61	25.36	113.22	505.52	771.84	11.94	18.23	81.40	363.43	554.90
19	16.48	25.16	112.35	501.59	765.84	11.85	18.09	80.77	360.61	550.59
20	16.36	24.97	111.49	497.78	760.02	11.76	17.95	80.15	357.87	546.40
21	16.23	24.78	110.66	494.05	754.32	11.67	17.82	79.55	355.19	542.30
22	16.11	24.60	109.82	490.31	748.61	11.58	17.68	78.95	352.49	538.20
23	15.99	24.41	108.98	486.58	742.93	11.49	17.55	78.35	349.82	534.11
24	15.87	24.22	108.16	482.89	737.29	11.41	17.42	77.76	347.16	530.06
25	15.74	24.03	107.31	479.11	731.51	11.32	17.28	77.15	344.44	525.90
26	15.61	23.84	106.43	475.16	725.49	11.22	17.14	76.51	341.61	521.58
27	15.48	23.63	105.50	471.04	719.19	11.13	16.99	75.85	338.64	517.05
28	15.33	23.41	104.51	466.63	712.46	11.02	16.83	75.14	335.47	512.21
29	15.18	23.17	103.46	461.95	705.31	10.91	16.66	74.38	332.11	507.07
30	15.02	22.94	102.41	457.22	698.09	10.80	16.49	73.62	328.71	501.88
31	14.88	22.71	101.40	452.74	691.26	10.69	16.33	72.90	325.49	496.96
32	14.74	22.51	100.48	448.63	684.98	10.60	16.18	72.24	322.53	492.45
33	14.62	22.33	99.69	445.07	679.55	10.51	16.05	71.67	319.98	488.54
34	14.53	22.18	99.04	442.17	675.12	10.44	15.95	71.20	317.89	485.36
35	14.45	22.07	98.53	439.90	671.65	10.39	15.87	70.83	316.26	482.87
36	14.40	21.98	98.15	438.21	669.07	10.35	15.80	70.56	315.04	481.01
37	14.35	21.92	97.85	436.88	667.03	10.32	15.76	70.35	314.08	479.55
38	14.31	21.85	97.58	435.65	665.16	10.29	15.71	70.15	313.20	478.21
39	14.28	21.80	97.31	434.47	663.36	10.26	15.67	69.96	312.35	476.91
40	14.24	21.74	97.07	433.39	661.71	10.24	15.63	69.79	311.58	475.72
41	14.21	21.69	96.84	432.37	660.16	10.21	15.59	69.62	310.85	474.61
42	14.17	21.64	96.61	431.32	658.55	10.19	15.56	69.45	310.09	473.45
43	14.13	21.57	96.32	430.02	656.57	10.16	15.51	69.24	309.16	472.03
44	14.07	21.48	95.88	428.10	653.63	10.11	15.44	68.93	307.77	469.91
45	13.98	21.34	95.29	425.45	649.58	10.05	15.34	68.51	305.87	467.00
46	13.87	21.18	94.55	422.15	644.54	9.97	15.22	67.98	303.49	463.38
47	13.74	20.98	93.67	418.20	638.52	9.88	15.08	67.34	300.66	459.05
48	13.59	20.76	92.67	413.74	631.71	9.77	14.92	66.62	297.45	454.16
49	13.44	20.52	91.61	409.02	624.50	9.66	14.75	65.86	294.06	448.97
50	13.29	20.28	90.57	404.36	617.38	9.55	14.58	65.11	290.70	443.85
51	13.15	20.07	89.61	400.10	610.88	9.45	14.43	64.43	287.64	439.18
52	13.03	19.90	88.84	396.65	605.61	9.37	14.31	63.87	285.16	435.39
53	12.96	19.78	88.32	394.32	602.05	9.31	14.22	63.49	283.49	432.83
54	12.92	19.73	88.10	393.35	600.58	9.29	14.19	63.34	282.79	431.77
55	12.94	19.76	88.24	393.96	601.50	9.31	14.21	63.44	283.23	432.44
56	13.01	19.87	88.72	396.11	604.78	9.36	14.29	63.78	284.77	434.79
57	13.13	20.05	89.52	399.68	610.24	9.44	14.41	64.36	287.34	438.72
58	13.29	20.29	90.60	404.50	617.59	9.55	14.59	65.13	290.80	444.01
59	13.48	20.59	91.92	410.39	626.60	9.69	14.80	66.08	295.04	450.48
60	13.71	20.93	93.46	417.27	637.10	9.86	15.05	67.19	299.99	458.03
61	13.97	21.32	95.21	425.07	649.01	10.04	15.33	68.45	305.60	466.59
62	14.25	21.76	97.15	433.75	662.26	10.25	15.64	69.84	311.83	476.12
63	14.57	22.24	99.29	443.31	676.85	10.47	15.99	71.38	318.71	486.61
64	14.91	22.77	101.67	453.94	693.09	10.72	16.37	73.10	326.35	498.28
65	15.31	23.37	104.35	465.90	711.34	11.01	16.80	75.02	334.95	511.40
66	15.76	24.06	107.42	479.61	732.27	11.33	17.30	77.23	344.80	526.45
67	16.29	24.87	111.02	495.69	756.83	11.71	17.88	79.82	356.37	544.11
68	16.92	25.83	115.32	514.87	786.11	12.16	18.57	82.91	370.15	565.16
69	17.67	26.98	120.47	537.86	821.21	12.70	19.40	86.61	386.68	590.39
70	18.57	28.36	126.60	565.26	863.05	13.35	20.39	91.02	406.38	620.47
71	19.63	29.97	133.82	597.49	912.26	14.11	21.55	96.21	429.55	655.85
72	20.86	31.85	142.20	634.88	969.35	15.00	22.90	102.23	456.43	696.89
73	22.27	34.00	151.82	677.83	1034.92	16.01	24.45	109.15	487.31	744.03
74	23.88	36.46	162.80	726.88	1109.82	17.17	26.22	117.04	522.58	797.88
75	25.72	39.26	175.31	782.70	1195.05	18.49	28.23	126.03	562.71	859.15
76	27.80	42.44	189.49	846.03	1291.74	19.98	30.51	136.23	608.24	928.67
77	30.15	46.04	205.54	917.68	1401.13	21.68	33.10	147.77	659.75	1007.32
78	32.81	50.09	223.66	998.57	1524.64	23.59	36.01	160.79	717.90	1096.10
79	35.80	54.67	244.07	1089.71	1663.78	25.74	39.30	175.47	783.42	1196.14
80	39.17	59.80	267.00	1192.11	1820.14	28.16	42.99	191.96	857.04	1308.55
81	42.93	65.55	292.65	1306.62	1994.97	30.86	47.12	210.40	939.37	1434.24
82	47.13	71.96	321.29	1434.47	2190.17	33.88	51.73	230.98	1031.28	1574.58
83	51.82	79.12	353.24	1577.12	2407.98	37.25	56.88	253.95	1133.84	1731.17
84	57.04	87.10	388.86	1736.19	2650.84	41.01	62.62	279.57	1248.19	1905.77
85	62.87	95.99	428.56	1913.40	2921.42	45.20	69.01	308.10	1375.60	2100.29
86	69.35	105.88	472.74	2110.67	3222.61	49.86	76.12	339.87	1517.42	2316.83
87	76.56	116.89	521.87	2330.02	3557.52	55.04	84.03	375.19	1675.12	2557.60
88	84.55	129.10	576.40	2573.48	3929.24	60.79	92.81	414.39	1850.15	2824.85

Table A8. Percentiles of IP-10 levels estimated with the fitted GAMLSS model. Family distribution Power Exponential with age dependence modeled as a linear function for location and scale parameters (μ and σ respectively). Sex is added as a qualitative variable for location and scale parameters (μ and σ respectively).

Age (Years)	Percentile for Women					Percentile for Men
Age (Years)	2.5th	5th	50th	95th	97.5th	2.5th	5th	50th	95th	97.5th
18	5548.71	8410.15	36,226.93	156,048.41	236,521.53	5986.14	8352.22	26,900.80	86,642.06	120,888.05
19	5428.03	8252.57	35,934.03	156,467.01	237,886.30	5872.14	8213.35	26,683.31	86,687.94	121,250.29
20	5309.43	8097.29	35,643.50	156,899.26	239,283.49	5759.84	8076.29	26,467.57	86,739.39	121,623.60
21	5192.88	7944.30	35,355.31	157,345.35	240,713.75	5649.21	7940.99	26,253.57	86,796.47	122,008.15
22	5078.36	7793.56	35,069.46	157,805.49	242,177.76	5540.25	7807.46	26,041.31	86,859.21	122,404.13
23	4965.85	7645.06	34,785.91	158,279.90	243,676.22	5432.94	7675.67	25,830.76	86,927.68	122,811.70
24	4855.32	7498.78	34,504.66	158,768.79	245,209.86	5327.25	7545.61	25,621.91	87,001.94	123,231.05
25	4746.74	7354.68	34,225.69	159,272.37	246,779.40	5223.17	7417.26	25,414.75	87,082.03	123,662.38
26	4640.09	7212.75	33,948.96	159,790.89	248,385.62	5120.69	7290.60	25,209.27	87,168.02	124,105.87
27	4535.35	7072.97	33,674.48	160,324.56	250,029.29	5019.78	7165.63	25,005.45	87,259.97	124,561.72
28	4432.49	6935.31	33,402.22	160,873.64	251,711.21	4920.43	7042.32	24,803.27	87,357.94	125,030.14
29	4331.49	6799.74	33,132.15	161,438.37	253,432.22	4822.63	6920.66	24,602.74	87,462.00	125,511.33
30	4232.33	6666.26	32,864.27	162,019.01	255,193.17	4726.35	6800.63	24,403.82	87,572.21	126,005.52
31	4134.97	6534.83	32,598.56	162,615.80	256,994.93	4631.58	6682.22	24,206.51	87,688.63	126,512.92
32	4039.40	6405.43	32,334.99	163,229.03	258,838.41	4538.31	6565.42	24,010.79	87,811.34	127,033.77
33	3945.59	6278.04	32,073.56	163,858.97	260,724.54	4446.51	6450.20	23,816.66	87,940.41	127,568.29
34	3853.53	6152.64	31,814.24	164,505.89	262,654.28	4356.17	6336.56	23,624.10	88,075.92	128,116.72
35	3763.18	6029.21	31,557.01	165,170.10	264,628.61	4267.27	6224.47	23,433.09	88,217.93	128,679.32
36	3674.53	5907.72	31,301.87	165,851.88	266,648.54	4179.81	6113.93	23,243.63	88,366.52	129,256.33
37	3587.54	5788.16	31,048.79	166,551.54	268,715.13	4093.75	6004.91	23,055.70	88,521.78	129,848.02
38	3502.21	5670.50	30,797.75	167,269.40	270,829.44	4009.09	5897.41	22,869.29	88,683.79	130,454.65
39	3418.50	5554.73	30,548.75	168,005.77	272,992.60	3925.81	5791.41	22,684.39	88,852.63	131,076.50
40	3336.40	5440.81	30,301.75	168,761.00	275,205.73	3843.90	5686.89	22,500.98	89,028.38	131,713.84
41	3255.88	5328.73	30,056.76	169,535.42	277,470.02	3763.33	5583.84	22,319.06	89,211.14	132,366.98
42	3176.92	5218.47	29,813.74	170,329.38	279,786.68	3684.09	5482.24	22,138.61	89,401.00	133,036.21
43	3099.50	5110.01	29,572.69	171,143.25	282,156.95	3606.18	5382.09	21,959.61	89,598.04	133,721.83
44	3023.59	5003.33	29,333.59	171,977.38	284,582.14	3529.56	5283.36	21,782.06	89,802.37	134,424.16
45	2949.18	4898.41	29,096.43	172,832.16	287,063.57	3454.23	5186.05	21,605.95	90,014.09	135,143.52
46	2876.24	4795.22	28,861.18	173,707.99	289,602.59	3380.17	5090.13	21,431.26	90,233.29	135,880.25
47	2804.76	4693.74	28,627.83	174,605.26	292,200.64	3307.37	4995.60	21,257.99	90,460.08	136,634.69
48	2734.71	4593.97	28,396.37	175,524.38	294,859.16	3235.81	4902.43	21,086.11	90,694.56	137,407.19
49	2666.07	4495.87	28,166.78	176,465.78	297,579.65	3165.48	4810.63	20,915.63	90,936.85	138,198.11
50	2598.82	4399.43	27,939.04	177,429.90	300,363.66	3096.36	4720.17	20,746.52	91,187.07	139,007.82
51	2532.94	4304.62	27,713.15	178,417.17	303,212.78	3028.43	4631.03	20,578.78	91,445.31	139,836.71
52	2468.41	4211.44	27,489.08	179,428.07	306,128.68	2961.69	4543.21	20,412.40	91,711.71	140,685.17
53	2405.20	4119.85	27,266.83	180,463.07	309,113.03	2896.12	4456.70	20,247.36	91,986.38	141,553.60
54	2343.31	4029.83	27,046.37	181,522.65	312,167.60	2831.69	4371.47	20,083.66	92,269.46	142,442.43
55	2282.71	3941.38	26,827.70	182,607.31	315,294.21	2768.41	4287.52	19,921.27	92,561.06	143,352.07
56	2223.38	3854.47	26,610.79	183,717.58	318,494.70	2706.25	4204.83	19,760.21	92,861.32	144,282.97
57	2165.30	3769.08	26,395.64	184,853.97	321,771.02	2645.20	4123.38	19,600.44	93,170.38	145,235.58
58	2108.45	3685.19	26,182.22	186,017.03	325,125.16	2585.25	4043.18	19,441.97	93,488.38	146,210.37
59	2052.81	3602.79	25,970.53	187,207.33	328,559.16	2526.38	3964.19	19,284.78	93,815.46	147,207.81
60	1998.36	3521.85	25,760.56	188,425.44	332,075.16	2468.58	3886.42	19,128.86	94,151.76	148,228.41
61	1945.09	3442.36	25,552.28	189,671.96	335,675.35	2411.83	3809.84	18,974.20	94,497.44	149,272.67
62	1892.97	3364.30	25,345.68	190,947.48	339,361.99	2356.12	3734.44	18,820.79	94,852.65	150,341.11
63	1842.00	3287.64	25,140.76	192,252.65	343,137.42	2301.44	3660.22	18,668.62	95,217.55	151,434.27
64	1792.14	3212.38	24,937.49	193,588.11	347,004.06	2247.78	3587.15	18,517.68	95,592.30	152,552.70
65	1743.38	3138.49	24,735.87	194,954.52	350,964.40	2195.11	3515.23	18,367.96	95,977.08	153,696.98
66	1695.70	3065.96	24,535.87	196,352.57	355,021.04	2143.43	3444.45	18,219.45	96,372.05	154,867.70
67	1649.09	2994.77	24,337.50	197,782.97	359,176.64	2092.73	3374.78	18,072.14	96,777.40	156,065.45
68	1603.52	2924.89	24,140.72	199,246.44	363,433.96	2042.98	3306.22	17,926.03	97,193.30	157,290.87
69	1558.99	2856.32	23,945.54	200,743.72	367,795.85	1994.18	3238.76	17,781.09	97,619.94	158,544.59
70	1515.46	2789.04	23,751.94	202,275.59	372,265.27	1946.32	3172.38	17,637.33	98,057.51	159,827.27
71	1472.94	2723.02	23,559.90	203,842.84	376,845.27	1899.38	3107.07	17,494.73	98,506.21	161,139.61
72	1431.39	2658.26	23,369.41	205,446.28	381,539.00	1853.35	3042.82	17,353.28	98,966.25	162,482.29
73	1390.80	2594.73	23,180.46	207,086.76	386,349.74	1808.21	2979.61	17,212.97	99,437.83	163,856.05
74	1351.15	2532.42	22,993.05	208,765.15	391,280.86	1763.96	2917.45	17,073.80	99,921.16	165,261.62
75	1312.44	2471.30	22,807.14	210,482.33	396,335.87	1720.58	2856.30	16,935.76	100,416.47	166,699.79
76	1274.63	2411.38	22,622.74	212,239.23	401,518.38	1678.05	2796.17	16,798.83	100,923.98	168,171.33
77	1237.72	2352.61	22,439.83	214,036.79	406,832.15	1636.38	2737.04	16,663.01	101,443.93	169,677.08
78	1201.70	2295.01	22,258.40	215,875.99	412,281.06	1595.54	2678.89	16,528.28	101,976.55	171,217.86
79	1166.53	2238.53	22,078.44	217,757.84	417,869.12	1555.51	2621.72	16,394.65	102,522.09	172,794.55
80	1132.22	2183.17	21,899.93	219,683.37	423,600.51	1516.30	2565.52	16,262.10	103,080.79	174,408.05
81	1098.73	2128.92	21,722.87	221,653.66	429,479.54	1477.89	2510.27	16,130.61	103,652.93	176,059.28
82	1066.07	2075.75	21,547.23	223,669.81	435,510.67	1440.27	2455.96	16,000.19	104,238.77	177,749.19
83	1034.20	2023.66	21,373.02	225,732.96	441,698.55	1403.42	2402.58	15,870.83	104,838.58	179,478.77
84	1003.13	1972.62	21,200.21	227,844.27	448,047.96	1367.33	2350.12	15,742.51	105,452.64	181,249.04
85	972.82	1922.61	21,028.81	230,004.97	454,563.89	1331.99	2298.57	15,615.23	106,081.24	183,061.03
86	943.28	1873.64	20,858.78	232,216.30	461,251.49	1297.39	2247.92	15,488.98	106,724.69	184,915.84
87	914.48	1825.67	20,690.14	234,479.55	468,116.12	1263.52	2198.15	15,363.75	107,383.28	186,814.57
88	886.41	1778.69	20,522.85	236,796.04	475,163.32	1230.37	2149.26	15,239.53	108,057.34	188,758.38

Table A9. Percentiles of MCP-1 levels estimated with the fitted GAMLSS model. Family distribution Normal with age linear dependence for location parameter (μ).

Age (Years)	Percentile
Age (Years)	2.5th	5th	50th	95th	97.5th
18	13.65	21.22	211.86	2115.55	3287.57
19	13.97	21.71	216.75	2164.31	3363.35
20	14.29	22.21	221.74	2214.19	3440.86
21	14.62	22.72	226.85	2265.22	3520.17
22	14.96	23.24	232.08	2317.43	3601.30
23	15.30	23.78	237.43	2370.84	3684.30
24	15.65	24.33	242.90	2425.49	3769.22
25	16.01	24.89	248.50	2481.39	3856.09
26	16.38	25.46	254.23	2538.58	3944.96
27	16.76	26.05	260.09	2597.09	4035.89
28	17.15	26.65	266.08	2656.94	4128.90
29	17.54	27.26	272.21	2718.18	4224.07
30	17.95	27.89	278.49	2780.83	4321.42
31	18.36	28.53	284.91	2844.92	4421.02
32	18.78	29.19	291.47	2910.49	4522.92
33	19.22	29.86	298.19	2977.57	4627.16
34	19.66	30.55	305.06	3046.20	4733.80
35	20.11	31.25	312.09	3116.40	4842.91
36	20.58	31.98	319.29	3188.23	4954.53
37	21.05	32.71	326.65	3261.71	5068.72
38	21.54	33.47	334.18	3336.89	5185.54
39	22.03	34.24	341.88	3413.80	5305.05
40	22.54	35.03	349.76	3492.48	5427.32
41	23.06	35.83	357.82	3572.97	5552.41
42	23.59	36.66	366.06	3655.32	5680.38
43	24.13	37.50	374.50	3739.57	5811.30
44	24.69	38.37	383.13	3825.75	5945.24
45	25.26	39.25	391.96	3913.93	6082.26
46	25.84	40.16	401.00	4004.14	6222.45
47	26.44	41.08	410.24	4096.42	6365.86
48	27.05	42.03	419.69	4190.84	6512.58
49	27.67	43.00	429.37	4287.43	6662.68
50	28.31	43.99	439.26	4386.24	6816.24
51	28.96	45.00	449.39	4487.34	6973.34
52	29.63	46.04	459.75	4590.76	7134.06
53	30.31	47.10	470.34	4696.56	7298.48
54	31.01	48.19	481.18	4804.81	7466.70
55	31.72	49.30	492.27	4915.55	7638.79
56	32.45	50.44	503.62	5028.84	7814.84
57	33.20	51.60	515.22	5144.75	7994.96
58	33.97	52.79	527.10	5263.32	8179.23
59	34.75	54.00	539.25	5384.63	8367.74
60	35.55	55.25	551.68	5508.73	8560.60
61	36.37	56.52	564.39	5635.70	8757.90
62	37.21	57.82	577.40	5765.59	8959.75
63	38.07	59.16	590.71	5898.47	9166.25
64	38.94	60.52	604.32	6034.42	9377.51
65	39.84	61.92	618.25	6173.50	9593.64
66	40.76	63.34	632.50	6315.78	9814.76
67	41.70	64.80	647.08	6461.35	10,040.96
68	42.66	66.30	661.99	6610.27	10,272.39
69	43.64	67.82	677.25	6762.62	10,509.14
70	44.65	69.39	692.86	6918.48	10,751.35
71	45.68	70.99	708.83	7077.94	10,999.15
72	46.73	72.62	725.16	7241.07	11,252.65
73	47.81	74.30	741.88	7407.96	11,512.00
74	48.91	76.01	758.97	7578.70	11,777.33
75	50.04	77.76	776.47	7753.37	12,048.77
76	51.19	79.55	794.36	7932.07	12,326.47
77	52.37	81.39	812.67	8114.88	12,610.56
78	53.58	83.26	831.40	8301.91	12,901.21
79	54.81	85.18	850.56	8493.25	13,198.55
80	56.08	87.14	870.17	8689.00	13,502.75
81	57.37	89.15	890.22	8889.27	13,813.96
82	58.69	91.21	910.74	9094.14	14,132.34
83	60.04	93.31	931.73	9303.74	14,458.06
84	61.43	95.46	953.21	9518.17	14,791.29
85	62.84	97.66	975.17	9737.55	15,132.19
86	64.29	99.91	997.65	9961.98	15,480.95
87	65.77	102.21	1020.64	10,191.58	15,837.76
88	67.29	104.57	1044.17	10,426.47	16,202.78

Table A10. Percentiles of VEGF levels estimated with the fitted GAMLSS model. Family distribution Power Exponential with age dependence modeled as a cubic spline for location parameter and as a linear function for skewness and scale parameters (μ and ν respectively).

Age (Years)	Percentile
Age (Years)	2.5th	5th	50th	95th	97.5th
18	23.58	31.23	124.00	492.38	652.17
19	24.37	32.09	127.31	504.99	665.03
20	25.19	32.99	130.71	517.91	678.20
21	26.04	33.92	134.22	531.18	691.73
22	26.93	34.88	137.86	544.87	705.72
23	27.86	35.89	141.65	559.09	720.32
24	28.83	36.95	145.61	573.88	735.55
25	29.85	38.06	149.79	589.43	751.66
26	30.94	39.26	154.25	606.05	769.04
27	32.11	40.55	159.05	623.88	787.83
28	33.34	41.90	164.08	642.56	807.60
29	34.60	43.29	169.22	661.55	827.63
30	35.87	44.68	174.37	680.49	847.49
31	37.14	46.06	179.42	699.00	866.72
32	38.39	47.39	184.31	716.78	884.96
33	39.59	48.68	188.96	733.56	901.89
34	40.74	49.89	193.32	749.11	917.25
35	41.83	51.01	197.32	763.26	930.86
36	42.84	52.04	200.94	775.88	942.58
37	43.74	52.93	204.02	786.36	951.70
38	44.49	53.64	206.39	794.07	957.48
39	45.07	54.15	207.97	798.76	959.67
40	45.49	54.46	208.80	800.53	958.42
41	45.75	54.58	208.91	799.58	954.02
42	45.86	54.53	208.36	796.11	946.73
43	45.82	54.31	207.18	790.28	936.75
44	45.66	53.95	205.45	782.38	924.46
45	45.41	53.49	203.37	773.23	910.83
46	45.12	52.99	201.17	763.65	896.86
47	44.83	52.50	198.97	754.13	883.10
48	44.55	52.02	196.87	745.02	869.95
49	44.30	51.58	194.91	736.53	857.66
50	44.07	51.18	193.14	728.76	846.34
51	43.89	50.84	191.56	721.80	836.05
52	43.75	50.55	190.20	715.70	826.89
53	43.65	50.30	189.03	710.34	818.66
54	43.56	50.09	187.97	705.47	811.09
55	43.49	49.88	186.97	700.84	803.88
56	43.40	49.67	185.97	696.22	796.76
57	43.30	49.44	184.88	691.37	789.45
58	43.16	49.17	183.67	686.05	781.69
59	42.96	48.85	182.26	680.06	773.24
60	42.70	48.45	180.60	673.18	763.85
61	42.36	47.97	178.64	665.20	753.31
62	41.93	47.40	176.32	655.93	741.38
63	41.40	46.71	173.60	645.22	727.91
64	40.76	45.90	170.45	632.96	712.79
65	40.01	44.98	166.88	619.18	696.04
66	39.17	43.95	162.95	604.14	677.96
67	38.25	42.85	158.76	588.15	658.92
68	37.29	41.70	154.39	571.54	639.27
69	36.28	40.52	149.90	554.57	619.32
70	35.26	39.31	145.35	537.37	599.20
71	34.21	38.09	140.74	520.01	578.98
72	33.15	36.85	136.08	502.52	558.71
73	32.07	35.60	131.40	484.96	538.43
74	30.98	34.35	126.70	467.38	518.21
75	29.88	33.09	122.00	449.84	498.11
76	28.78	31.83	117.32	432.39	478.18
77	27.69	30.58	112.67	415.08	458.48
78	26.60	29.34	108.06	397.97	439.05
79	25.52	28.12	103.51	381.08	419.93
80	24.45	26.91	99.03	364.47	401.18
81	23.39	25.72	94.63	348.17	382.82
82	22.36	24.56	90.33	332.26	364.95
83	21.35	23.43	86.15	316.82	347.63
84	20.38	22.34	82.12	301.90	330.94
85	19.43	21.28	78.23	287.56	314.93
86	18.53	20.27	74.51	273.84	299.63
87	17.66	19.31	70.96	260.76	285.06
88	16.84	18.40	67.59	248.32	271.23

Appendix C. Reference Intervals Classification

This appendix shows the complete classification of all individuals of validation and utility samples based on the cytokine reference limits estimated by the fitted generalized additive models for location, scale and shape (GAMLSS) models.

Table A11. Complete classification of validation sample by reference intervals (RIs) estimated with the fitted GAMLSS model.

Age	Sex	OSDI	T-BUT	Schirmer	Corneal Stain.	EGF		Fractalkine		IL-1RA		IL-8		IP-10		MCP-1		VEGF
Age	Sex	OSDI	T-BUT	Schirmer	Corneal Stain.	Level	Class	Level	Class	Level	Class	Level	Class	Level	Class	Level	Class	Level	Class
27	M	2.1	10.0	30.0	0	1940	Normal	424	Normal	2020	Normal	293	Normal	54,000	Normal	NA	NA	555	Normal
29	F	4.5	9.7	35.0	0	2660	Normal	630	Normal	2700	Normal	429	Normal	48,600	Normal	NA	NA	58	Normal
25	F	4.2	8.3	12.0	0	2110	Normal	237	Normal	2090	Normal	186	Normal	62,800	Normal	NA	NA	58	Normal
29	M	0.0	10.0	8.0	0	1770	Normal	890	Normal	2330	Normal	155	Normal	29,200	Normal	NA	NA	58	Normal
56	M	0.0	10.3	5.0	0	419	Normal	NA	NA	8400	High	122.5	Normal	67,500	Normal	NA	NA	58	Normal
24	M	0.0	15.0	9.0	0	880	Normal	469.5	Normal	3175	Normal	486	High	50,500	Normal	NA	NA	520	Normal
20	F	0.0	10.0	20.0	0	1890	Normal	221.5	Normal	5350	Normal	139	Normal	22,250	Normal	NA	NA	104	Normal
19	F	4.5	11.0	14.0	0	1110	Normal	560	Normal	NA	NA	98	Normal	31,450	Normal	NA	NA	429.5	Normal
20	M	4.2	9.0	15.0	0	1650	Normal	313	Normal	1935	Normal	133	Normal	38,350	Normal	NA	NA	104	Normal
18	M	4.2	16.0	22.0	0	1675	Normal	640	Normal	4620	Normal	421	High	54,500	Normal	NA	NA	460	Normal
19	F	10.4	18.0	14.0	0	1635	Normal	448.5	Normal	7950	High	465	Normal	41,950	Normal	NA	NA	465	Normal
20	F	0.0	13.0	6.0	0	449	Normal	510	Normal	590	Normal	190.5	Normal	41,850	Normal	NA	NA	715	High
30	F	4.2	15.0	30.0	0	520	Normal	685	Normal	258.5	Normal	385	Normal	51,500	Normal	NA	NA	364.5	Normal
20	F	8.3	18.0	35.0	0	865	Normal	235.5	Normal	4380	Normal	148	Normal	47,150	Normal	NA	NA	104	Normal
20	F	2.1	18.0	8.0	0	1335	Normal	560	Normal	275.5	Normal	340.5	Normal	23,850	Normal	NA	NA	354.5	Normal
20	F	8.3	8.0	6.0	0	235	Normal	416.5	Normal	275.5	Normal	69.5	Normal	13,700	Normal	NA	NA	333.5	Normal
20	F	8.3	25.0	15.0	0	2565	Normal	890	High	4395	Normal	397	Normal	52,500	Normal	NA	NA	104	Normal
20	F	8.3	15.0	10.0	0	1105	Normal	349	Normal	155	Normal	92.5	Normal	39,700	Normal	NA	NA	354.5	Normal
24	M	2.1	11.0	35.0	0	870	Normal	427.5	Normal	1600	Normal	116	Normal	30,800	Normal	NA	NA	104	Normal
24	F	4.2	11.0	35.0	0	805	Normal	325	Normal	4600	Normal	429.5	Normal	63,500	Normal	NA	NA	645	High
47	F	2.1	13.0	33.0	0	170	Normal	368	Normal	296	Normal	30	Normal	4790	Normal	NA	NA	463	Normal
66	M	0.0	13.0	12.0	0	1190	Normal	1040	Normal	1120	Normal	145	Normal	7760	Normal	NA	NA	550	Normal
53	M	10.4	9.0	6.0	0	1560	Normal	1490	High	849	Normal	80	Normal	19,200	Normal	NA	NA	242	Normal
62	M	4.2	9.0	10.0	0	406	Normal	1050	Normal	3320	Normal	150	Normal	59,400	Normal	NA	NA	489	Normal
25	F	0.0	12.0	20.0	0	1280	Normal	187	Normal	545	Normal	296	Normal	9060	Normal	NA	NA	108	Normal
62	M	0.0	12.0	10.0	0	566	Normal	187	Normal	2420	Normal	235	Normal	3930	Normal	NA	NA	108	Normal
18	F	0.0	12.0	25.0	0	597	Normal	770	Normal	1190	Normal	357	Normal	40,400	Normal	NA	NA	108	Normal
21	M	0.0	10.0	20.0	0	275	Normal	187	Normal	115	Normal	102	Normal	7910	Low	NA	NA	108	Normal
56	M	0.0	12.0	20.0	0	408	Normal	187	Normal	186	Normal	262	Normal	16,900	Normal	NA	NA	108	Normal
62	M	0.0	10.0	15.0	0	518	Normal	679	Normal	5890	High	306	Normal	18,600	Normal	NA	NA	108	Normal
58	F	0.0	12.0	13.0	0	164	Normal	418	Normal	2620	Normal	199	Normal	24,600	Normal	NA	NA	108	Normal
54	M	0.0	11.7	22.0	0	544	Normal	261	Normal	5525	Normal	240	Normal	23,300	Normal	NA	NA	630	Normal
31	M	9.1	16.0	35.0	0	171	Normal	1050	Normal	4780	Normal	180	Normal	NA	NA	524	Normal	NA	NA
24	F	11.4	11.7	35.0	0	1670	Normal	NA	NA	2570	Normal	253	Normal	27,000	Normal	1200	Normal	NA	NA
19	M	8.3	8.7	30.0	0	383	Normal	NA	NA	740	Normal	29.3	Normal	23,500	Normal	143	Normal	NA	NA
24	F	6.3	7.3	35.0	0	1670	Normal	120	Normal	1120	Normal	92.2	Normal	NA	NA	126	Normal	NA	NA
36	M	5.6	11.3	19.0	0	864	Normal	318	Normal	5700	Normal	95.6	Normal	NA	NA	1500	Normal	NA	NA
21	F	11.1	8.0	5.0	0	1530	Normal	NA	NA	1610	Normal	256	Normal	23,400	Normal	189	Normal	NA	NA
20	F	7.5	8.0	22.0	0	551	Normal	225	Normal	4510	Normal	145	Normal	NA	NA	141	Normal	NA	NA
27	M	9.1	7.3	30.0	0	999	Normal	800	Normal	3660	Normal	101	Normal	NA	NA	166	Normal	NA	NA

M = Male; F = Female; OSDI = Ocular surface index score; T-BUT = Fluorescein tear break-up time; NA = Not Available.

Table A12. Complete classification of utility sample by reference intervals (RIs) estimated with the fitted GAMLSS model.

Age	Sex	OSDI	T-BUT	Schirmer	Corneal Stain.	EGF		Fractalkine		IL-1RA		IL-8		IP-10		MCP-1		VEGF
Age	Sex	OSDI	T-BUT	Schirmer	Corneal Stain.	Level	Class	Level	Class	Level	Class	Level	Class	Level	Class	Level	Class	Level	Class
59	F	90.0	2.3	0.0	4.0	1040	Normal	1220	High	96,500	High	1780	High	110,000	Normal	NA	NA	1490	High
55	M	45.0	1.0	0.0	4.0	937	Normal	NA	NA	100,000	High	16,900	High	3300	Low	NA	NA	2100	High
47	F	66.7	0.7	0.0	3.0	117	Normal	615	Normal	45,400	High	327	Normal	244	Low	NA	NA	219	Normal
69	M	57.5	2.3	2.0	3.0	104	Normal	1050	Normal	48,300	High	1480	High	3230	Low	NA	NA	862	High
55	M	75.0	1.0	2.0	3.0	175	Normal	1000	Normal	43,100	High	2390	High	2430	Low	NA	NA	1050	High
54	F	75.0	1.0	1.0	3.0	292	Normal	2660	High	47,900	High	13,400	High	38,600	Normal	NA	NA	3050	High
60	F	78.1	1.0	0.0	3.0	468	Normal	1750	High	65,300	High	1370	High	40,400	Normal	NA	NA	699	High
58	F	45.0	2.0	2.0	3.0	27	Low	3890	High	88,800	High	3760	High	149,000	Normal	NA	NA	805	High
72	F	85.0	2.0	2.0	4.0	303	Normal	860	High	7420	High	1640	High	3800	Normal	NA	NA	1020	High
53	M	85.0	1.3	1.0	3.0	1040	Normal	1290	High	9390	High	19,700	High	12	Low	NA	NA	58	Normal
75	F	34.4	1.0	0.0	3.0	360	Normal	1490	High	10,400	High	4830	High	29,300	Normal	NA	NA	254	Normal
51	F	50.0	0.7	1.0	3.0	27	Low	5950	High	69,500	High	4470	High	32,500	Normal	NA	NA	2870	High
49	F	61.4	1.0	1.0	4.0	37.3	Low	NA	NA	1160	Normal	582	High	9250	Normal	12.3	Low	NA	NA

M = Male; F = Female; OSDI = Ocular surface index score; T-BUT = Fluorescein tear break-up time; NA = Not Available.

References

Monastero, R.N.; Pentyala, S. Cytokines as biomarkers and their respective clinical cutoff levels. Int. J. Inflamm. 2017, 4309485. [Google Scholar] [CrossRef] [PubMed] [Green Version]
von Thun Und Hohenstein-Blaul, N.; Funke, S.; Grus, F.H. Tears as a source of biomarkers for ocular and systemic diseases. Exp. Eye Res. 2013, 117, 126–137. [Google Scholar] [CrossRef] [PubMed]
Hagan, S.; Martin, E.; Enríquez-de-Salamanca, A. Tear fluid biomarkers in ocular and systemic disease: Potential use for predictive, preventive and personalised medicine. EPMA J. 2016, 7, 15. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Tamhane, M.; Cabrera-Ghayouri, S.; Abelian, G.; Viswanath, V. Review of biomarkers in ocular matrices: Challenges and opportunities. Pharm. Res. 2019, 36, 40. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Roy, N.S.; Wei, Y.; Kuklinski, E.; Asbell, P.A. The growing need for validated biomarkers and endpoints for dry eye clinical research. Investig. Ophthalmol. Vis. Sci. 2017, 58, BIO1–BIO19. [Google Scholar] [CrossRef] [Green Version]
Soria, J.; Durán, J.A.; Etxebarria, J.; Merayo, J.; González, N.; Reigada, R.; García, I.; Acera, A.; Suárez, T. Tear proteome and protein network analyses reveal a novel pentamarker panel for tear film characterization in dry eye and meibomian gland dysfunction. J. Proteom. 2013, 78, 94–112. [Google Scholar] [CrossRef] [PubMed]
Versura, P.; Bavelloni, A.; Grillini, M.; Fresina, M.; Campos, E.C. Diagnostic performance of a tear protein panel in early dry eye. Mol. Vis. 2013, 19, 1247–1257. [Google Scholar]
Pinto-Fraga, J.; Enríquez-de-Salamanca, A.; Calonge, M.; González-García, M.J.; López-Miguel, A.; López-de la Rosa, A.; García-Vázquez, C.; Calder, V.; Stern, M.E.; Fernández, I. Severity, therapeutic, and activity tear biomarkers in dry eye disease: An analysis from a phase III clinical trial. Ocul. Surf. 2018, 6, 368–376. [Google Scholar] [CrossRef]
Stern, M.E.; Schaumburg, C.S.; Pflugfelder, S.C. Dry eye as a mucosal autoimmune disease. Int. Rev. Immunol. 2013, 32, 19–41. [Google Scholar] [CrossRef]
Craig, J.P.; Nichols, K.K.; Akpek, E.K.; Caffery, B.; Dua, H.S.; Joo, C.K.; Liu, Z.; Nelson, J.D.; Nichols, J.J.; Tsubota, K.; et al. TFOS DEWS II Definition and Classification Report. Ocul. Surf. 2017, 15, 276–283. [Google Scholar] [CrossRef]
Zhou, L.; Beuerman, R.W. Tear analysis in ocular surface diseases. Prog. Retin. Eye Res. 2012, 31, 527–550. [Google Scholar] [CrossRef]
Micera, A.; Zazzo, A.D.; Esposito, G.; Longo, R.; Foulsham, W.; Sacco, R.; Sgrulletta, R.; Bonini, S. Age-related changes to human tear composition. Investig. Ophthalmol. Vis. Sci. 2018, 59, 2024–2031. [Google Scholar] [CrossRef] [Green Version]
Carreño, E.; Enríquez-de-Salamanca, A.; Tesón, M.; García-Vázquez, C.; Stern, M.E.; Whitcup, S.M.; Calonge, M. Cytokine and chemokine levels in tears from healthy subjects. Acta Ophthalmol. 2010, 88, e250–e258. [Google Scholar] [CrossRef]
LaFrance, M.W.; Kehinde, L.E.; Fullard, R.J. Multiple cytokine analysis in human tears: An optimized procedure for cytometric bead-based assay. Curr. Eye Res. 2008, 33, 525–544. [Google Scholar] [CrossRef]
McGill, J.I.; Liakos, G.M.; Goulding, N.; Seal, D.V. Normal tear protein profiles and age-related changes. Br. J. Ophthalmol. 1984, 68, 316–320. [Google Scholar] [CrossRef] [Green Version]
Nättinen, J.; Jylhä, A.; Aapola, U.; Mäkinen, P.; Beuerman, R.; Pietilä, J.; Vaajanen, A.; Uusitalo, H. Age-associated changes in human tear proteome. Clin. Proteom. 2019, 16, 11. [Google Scholar] [CrossRef] [PubMed]
Rigby, R.A.; Stasinopoulos, D.M. Generalized additive models for location, scale and shape (with discussion). Appl. Stat. 2005, 54, 507–554. [Google Scholar] [CrossRef] [Green Version]
Borghi, E.; de Onis, M.; Garza, C.; Van den Broeck, J.; Frongillo, E.A.; Grummer-Strawn, L.; Van Buuren, S.; Pan, H.; Molinari, L.; Martorell, R.; et al. Construction of the World Health Organization child growth standards: Selection of methods for attained growth curves. Stat. Med. 2006, 25, 247–265. [Google Scholar] [CrossRef] [PubMed]
Chamitava, L.; Garcia-Larsen, V.; Cazzoletti, L.; Degan, P.; Pasini, A.; Bellisario, V.; Corsico, A.G.; Nicolis, M.; Olivieri, M.; Pirina, P.; et al. Determination of adjusted reference intervals of urinary biomarkers of oxidative stress in healthy adults using GAMLSS models. PLoS ONE 2018, 13, e0206176. [Google Scholar] [CrossRef] [PubMed]
Flatley, C.; Kumar, S.; Greer, R.M. Reference centiles for the middle cerebral artery and umbilical artery pulsatility index and cerebro-placental ratio from a low-risk population—A Generalised Additive Model for Location, Shape and Scale (GAMLSS) approach. J. Matern.-Fetal Neonatal. Med. 2019, 32, 2338–2345. [Google Scholar] [CrossRef] [PubMed]
Enríquez-de-Salamanca, A.; Castellanos, E.; Stern, M.E.; Fernández, I.; Carreño, E.; García-Vázquez, C.; Herreras, J.M.; Calonge, M. Tear cytokine and chemokine analysis and clinical correlations in evaporative-type dry eye disease. Mol. Vis. 2010, 16, 862–873. [Google Scholar]
Carreño, E.; Portero, A.; Herreras, J.M.; García-Vázquez, C.; Whitcup, S.M.; Stern, M.E.; Calonge, M.; Enríquez-de-Salamanca, A. Cytokine and chemokine tear levels in patients with uveitis. Acta Ophthalmol. 2017, 95, e405–e414. [Google Scholar] [CrossRef] [Green Version]
Helsel, D.R.; Cohn, T.A. Estimation of descriptive statistics for multiply censored water quality data. Water Resour. Res. 1988, 24, 1997–2004. [Google Scholar] [CrossRef]
Akaike, H. Information measures and model selection. Int. Stat. Inst. 1983, 50, 277–291. [Google Scholar]
Dunn, P.K.; Smyth, G.K. Randomised quantile residuals. J. Comput. Graph. Stat. 1996, 5, 236–244. [Google Scholar] [CrossRef]
van Buuren, S.; Fredriks, M. Worm plot: A simple diagnostic device for modelling growth reference curves. Stat. Med. 2001, 20, 1259–1277. [Google Scholar] [CrossRef]
Wolffsohn, J.S.; Arita, R.; Chalmers, R.; Djalilian, A.; Dogru, M.; Dumbleton, K.; Gupta, P.K.; Karpecki, P.; Lazreg, S.; Pult, H.; et al. TFOS Dry Eye workshop (DEWS) II Diagnostic methodology report. Ocul. Surf. 2017, 15, 539–574. [Google Scholar] [CrossRef]
CLSI. Defining, Establishing, and Verifying Reference Intervals in the Clinical Laboratory; Approved Guideline, 3rd ed.; CLSI Document EP28-A3c; Clinical and Laboratory Standards Institute: Wayne, PA, USA, 2008. [Google Scholar]
Zhou, X.; Fragala, M.S.; McElhaney, J.E.; Kuchel, G.A. Conceptual and methodological issues relevant to cytokine and inflammatory marker measurements in clinical research. Curr. Opin. Clin. Nutr. Metab. Care 2010, 13, 541–547. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Genser, B.; Cooper, P.J.; Yazdanbakhsh, M.; Barreto, M.L.; Rodrigues, L.C. A guide to modern statistical analysis of immunological data. BMC Immunol. 2007, 8, 27. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Rea, I.M.; Gibson, D.S.; McGilligan, V.; McNerlan, S.E.; Alexander, H.D.; Ross, O.A. Age and age-related diseases: Role of inflammation triggers and cytokines. Front. Immunol. 2018, 9, 586. [Google Scholar] [CrossRef] [PubMed]
Franceschi, C.; Bonafè, M.; Valensin, S.; Olivieri, F.; De Luca, M.; Ottaviani, E.; De Benedictis, G. Inflamm-aging. An evolutionary perspective on immunosenescence. Ann. N. Y. Acad. Sci. 2000, 908, 244–254. [Google Scholar] [CrossRef]
Xia, S.; Zhang, X.; Zheng, S.; Khanabdali, R.; Kalionis, B.; Wu, J.; Wan, W.; Tai, X. An update on inflamm-aging: Mechanisms, prevention, and treatment. J. Immunol. Res. 2016, 2016, 8426874. [Google Scholar] [CrossRef]
Krabbe, K.S.; Pedersen, M.; Bruunsgaard, H. Inflammatory mediators in the elderly. Exp. Gerontol. 2004, 39, 687–699. [Google Scholar] [CrossRef]
Antonelli, A.; Rotondi, M.; Fallahi, P.; Ferrari, S.M.; Paolicchi, A.; Romagnani, P.; Serio, M.; Ferrannini, E. Increase of CXC chemokine CXCL10 and CC chemokine CCL2 serum levels in normal ageing. Cytokine 2006, 34, 32–38. [Google Scholar] [CrossRef]
Hu, W.T.; Howell, J.C.; Ozturk, T.; Gangishetti, U.; Kollhoff, A.L.; Hatcher-Martin, J.M.; Anderson, A.M.; Tyor, W.R. CSF cytokines in aging, multiple sclerosis, and dementia. Front. Immunol. 2019, 10, 480. [Google Scholar] [CrossRef] [Green Version]
Valiathan, R.; Ashman, M.; Asthana, D. Effects of ageing on the immune system: Infants to elderly. Scand. J. Immunol. 2016, 83, 255–266. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Wieczorowska-Tobis, K.; Niemir, Z.I.; Podkówka, R.; Korybalska, K.; Mossakowska, M.; Breborowicz, A. Can an increased level of circulating IL-8 be a predictor of human longevity? Med. Sci. Monit. 2006, 12, CR118-21. [Google Scholar] [PubMed]
Cavallone, L.; Bonafè, M.; Olivieri, F.; Cardelli, M.; Marchegiani, F.; Giovagnetti, S.; Di Stasio, G.; Giampieri, C.; Mugianesi, E.; Stecconi, R.; et al. The role of IL-1 gene cluster in longevity: A study in Italian population. Mech. Ageing Dev. 2003, 124, 533–538. [Google Scholar] [CrossRef]
Jylhä, M.; Paavilainen, P.; Lehtimäki, T.; Goebeler, S.; Karhunen, P.J.; Hervonen, A.; Hurme, M. Interleukin-1 receptor antagonist, interleukin-6, and C-reactive protein as predictors of mortality in nonagenarians: The vitality 90+ study. J. Gerontol. A Biol. Sci. Med. Sci. 2007, 62, 1016–1021. [Google Scholar] [CrossRef] [Green Version]
Nguyen, V.; Mendelsohn, A.; Larrick, J.W. Interleukin-7 and immunosenescence. J. Immunol. Res. 2017, 2017, 4807853. [Google Scholar] [CrossRef] [Green Version]
Nasi, M.; Troiano, L.; Lugli, E.; Pinti, M.; Ferraresi, R.; Monterastelli, E.; Mussi, C.; Salvioli, G.; Franceschi, C.; Cossarizza, A. Thymic output and functionality of the IL-7/IL-7 receptor system in centenarians: Implications for the neolymphogenesis at the limit of human life. Aging Cell 2006, 5, 167–175. [Google Scholar] [CrossRef] [Green Version]
Ruggiero, D.; Dalmasso, C.; Nutile, T.; Sorice, R.; Dionisi, L.; Aversano, M.; Bröet, P.; Leutenegger, A.L.; Bourgain, C.; Ciullo, M. Genetics of VEGF serum variation in human isolated populations of cilento: Importance of VEGF polymorphisms. PLoS ONE 2011, 6, e16982. [Google Scholar] [CrossRef]
Saldías, M.P.; Fernández, C.; Morgan, A.; Díaz, C.; Morales, D.; Jaña, F.; Gómez, A.; Silva, A.; Briceño, F.; Oyarzún, A.; et al. Aged blood factors decrease cellular responses associated with delayed gingival wound repair. PLoS ONE 2017, 12, e0184189. [Google Scholar] [CrossRef] [Green Version]
Scheubel, R.J.; Zorn, H.; Silber, R.E.; Kuss, O.; Morawietz, H.; Holtz, J.; Simm, A. Age-dependent depression in circulating endothelial progenitor cells in patients undergoing coronary artery bypass grafting. J. Am. Coll. Cardiol. 2003, 42, 2073–2080. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Chakraborty, A.; Chatterjee, M.; Twaalfhoven, H.; Del Campo Milan, M.; Teunissen, C.E.; Scheltens, P.; Fontijn, R.U.; van Der Flier, W.M.; de Vries, H.E. Vascular Endothelial Growth Factor remains unchanged in cerebrospinal fluid of patients with Alzheimer’s disease and vascular dementia. Alzheimers Res. Ther. 2018, 10, 58. [Google Scholar] [CrossRef] [PubMed]
Rübenhagen, R.; Schüttrumpf, J.P.; Stürmer, K.M.; Frosch, K.H. Interleukin-7 levels in synovial fluid increase with age and MMP-1 levels decrease with progression of osteoarthritis. Acta Orthop. 2012, 83, 59–64. [Google Scholar] [CrossRef] [Green Version]
Di Zazzo, A.; Micera, A.; Coassin, M.; Varacalli, G.; Foulsham, W.; De Piano, M.; Bonini, S. InflammAging at ocular surface: Clinical and biomolecular analyses in healthy volunteers. Investig. Ophthalmol. Vis. Sci. 2019, 60, 1769–1775. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Gallenga, C.E.; Parmeggiani, F.; Costagliola, C.; Sebastiani, A.; Gallenga, P.E. Inflammaging: Should this term be suitable for age related macular degeneration too? Inflamm. Res. 2014, 63, 105–107. [Google Scholar] [CrossRef] [PubMed]
Vohra, R.; Tsai, J.C.; Kolko, M. The role of inflammation in the pathogenesis of glaucoma. Surv. Ophthalmol. 2013, 58, 311–320. [Google Scholar] [CrossRef] [PubMed]
McMonnies, C.W. Glaucoma history and risk factors. J. Optom. 2017, 10, 71–78. [Google Scholar] [CrossRef] [Green Version]
Milan-Mattos, J.C.; Anibal, F.F.; Perseguini, N.M.; Minatel, V.; Rehder-Santos, P.; Castro, C.A.; Vasilceac, F.A.; Mattiello, S.M.; Faccioli, L.H.; Catai, A.M. Effects of natural aging and gender on pro-inflammatory markers. Braz. J. Med. Biol. Res. 2019, 52, e8392. [Google Scholar] [CrossRef] [Green Version]
Enríquez-de-Salamanca, A.; Calonge, M. Cytokines and chemokines in immune-based ocular surface inflammation. Expert Rev. Clin. Immunol. 2008, 4, 457–467. [Google Scholar] [CrossRef] [PubMed]
Enríquez-de-Salamanca, A.; Bonini, S.; Calonge, M. Molecular and cellular biomarkers in dry eye disease and ocular allergy. Curr. Opin. Allergy Clin. Immunol. 2012, 12, 523–533. [Google Scholar] [CrossRef] [PubMed]
Lam, H.; Bleiden, L.; de Paiva, C.S.; Farley, W.; Stern, M.E.; Pflugfelder, S.C. Tear cytokine profiles in dysfunctional tear syndrome. Am. J. Ophthalmol. 2009, 147, 198–205. [Google Scholar] [CrossRef] [PubMed] [Green Version]
Massingale, M.L.; Li, X.; Vallabhajosyula, M.; Chen, D.; Wei, Y.; Asbell, P.A. Analysis of inflammatory cytokines in the tears of dry eye patients. Cornea 2009, 28, 1023–1027. [Google Scholar] [CrossRef] [PubMed]
Choi, W.; Li, Z.; Oh, H.J.; Im, S.K.; Lee, S.H.; Park, S.H.; You, I.C.; Yoon, K.C. Expression of CCR5 and its ligands CCL3, −4, and −5 in the tear film and ocular surface of patients with dry eye disease. Curr. Eye Res. 2012, 37, 12–17. [Google Scholar] [CrossRef]
Huang, J.F.; Zhang, Y.; Rittenhouse, K.D.; Pickering, E.H.; McDowell, M.T. Evaluations of tear protein markers in dry eye disease: Repeatability of measurement and correlation with disease. Investig. Ophthalmol. Vis. Sci. 2012, 53, 4556–4564. [Google Scholar] [CrossRef] [Green Version]
Na, K.S.; Mok, J.W.; Kim, J.Y.; Rho, C.R.; Joo, C.K. Correlations between tear cytokines, chemokines, and soluble receptors and clinical severity of dry eye disease. Investig. Ophthalmol. Vis. Sci. 2012, 53, 5443–5450. [Google Scholar] [CrossRef] [Green Version]
Pflugfelder, S.C.; Jones, D.; Ji, Z.; Afonso, A. Monroy, D. Altered cytokine balance in the tear fluid and conjunctiva of patients with Sjögren’s syndrome keratoconjunctivitis sicca. Curr. Eye Res. 1999, 19, 201–211. [Google Scholar] [CrossRef]
Ohashi, Y.; Ishida, R.; Kojima, T.; Goto, E.; Matsumoto, Y.; Watanabe, K.; Ishida, N.; Nakata, K.; Takeuchi, T.; Tsubota, K. Abnormal protein profiles in tears with dry eye syndrome. Am. J. Ophthalmol. 2003, 136, 291–299. [Google Scholar] [CrossRef]
Cocho, L.; Fernández, I.; Calonge, M.; Sainz de la Maza, M.; Rovira, M.; Stern, M.E.; García-Vázquez, C.; Enríquez-de-Salamanca, A. Prehematopoietic stem cell transplantation tear cytokines as potential susceptibility biomarkers for ocular chronic Graft-Versus-Host disease. Investig. Ophthalmol. Vis. Sci. 2017, 58, 4836–4846. [Google Scholar] [CrossRef] [Green Version]
Gurumurthy, S.; Iyer, G.; Srinivasan, B.; Agarwal, S.; Angayarkanni, N. Ocular surface cytokine profile in chronic Stevens-Johnson syndrome and its response to mucous membrane grafting for lid margin keratinisation. Br. J. Ophthalmol. 2018, 102, 169–176. [Google Scholar] [CrossRef] [Green Version]
Hernandez-Molina, G.; Ruiz-Quintero, N.; Lima, G.; Hernández-Ramírez, D.; Llorente-Chavez, A.; Saavedra-González, V.; Jimenez-Soto, R.; Llorente, L. Dry eye in Sjögren’s Syndrome: Chemokine and Cytokine tear spectrum. Ann. Rheum. Dis. 2019, 78, 389. [Google Scholar] [CrossRef] [Green Version]
Than, T.P.; Guyette, N.; Tran, K.; Elder, K.; Sims, J.; Kehinde, L.; Shin, P.; Williams, L.; Bradley, J.; Fullard, R. Correlations between tear IP-10 and other biomarkers in normal and dry eye patients. Investig. Ophthalmol. Vis. Sci. 2011, 52, 3807. [Google Scholar]
Yoon, K.C.; Park, C.S.; You, I.C.; Choi, H.J.; Lee, K.H.; Im, S.K.; Park, H.Y.; Pflugfelder, S.C. Expression of CXCL9, -10, -11, and CXCR3 in the tear film and ocular surface of patients with dry eye syndrome. Investig. Ophthalmol. Vis. Sci. 2010, 51, 643–650. [Google Scholar] [CrossRef] [PubMed]

Figure 1. Stepwise algorithm for selecting relevant explanatory variables given a distribution family. Forward selection stage starts with the simplest model, the model with all the parameters fitted as constants, and tests the addition of sex and age, as a constant, a linear function or a cubic spline, using the Akaike Information Criterion (AIC). Backward elimination stage starts with the best model fitted in the forward stage, and tests the deletion of each explanatory variable using the AIC. The distributional parameters are: μ location; σ scale; ν skewness; and τ kurtosis.

Figure 2. Age- and sex-adjusted 90% reference interval (RI) for eight tear cytokine levels estimated with GAMLSS. Three estimated centile curves (5th, 50th, and 95th percentiles) of each tear cytokine level is represented. Observed cytokine levels are plotted over the range of ages. Circles and triangles represent female and male values respectively. Sex- dependent cytokines levels (Fractalkine, IL-7, IL-8, and IP-10) are shown as separate graph for females and males. Pg/mL = picograms per milliliter.

Table 1. Minimum detectable concentration and detection rate of the 18 cytokines analyzed in 157 tear samples.

	Minimum Detectable Concentration (pg/mL)	Detection Rate
	Minimum Detectable Concentration (pg/mL)	%	95% CI for %
EGF	3.2	100	97.02; 100
Fractalkine	3.83	87.9	81.51; 92.37
IFN-γ	0.1	1.3	0.22; 5
IL-10	0.3	36.3	28.9; 44.4
IL-12p70	0.4	5.1	2.39; 10.13
IL-15	0.4	29.9	23.03; 37.84
IL-17A	0.2	0	0; 2.98
IL-1β	0.4	19.1	13.45; 26.31
IL-1RA	3.2	100	97.02; 100
IL-2	0.3	0.6	0.03; 4.03
IL-6	0.3	21	15.1; 28.39
IL-7	1.08	87.3	80.78; 91.86
IL-8	0.2	96.8	92.34; 98.82
IP-10	3.2	100	97.02; 100
MCP-1	3.2	100	97.02; 100
TNF-α	0.1	1.9	0.49; 5.92
VEGF	3.73	82.2	75.08; 87.63
IL-23 (*)	1.01	10.2	6.02; 16.55

CI = Confidence interval; EGF = Epidermal Growth Factor; IFN = Interferon; IL = Interleukin; RA = Receptor Antagonist; IP = Induced Protein; MCP = Monocyte Chemotactic Protein; TNF = Tumor Necrosis Factor; VEGF = Vascular Endothelial Growth Factor. (*) Out of 157 subjects, 10 did not have data on IL-23.

Table 2. Description of tear cytokine levels by sex and age groups.

Age Range (n)	EGF	Fractalkine	IL-1RA	IL-7	IL-8	IP-10	MCP-1	VEGF
Age Range (n)	(pg/mL)	(pg/mL)	(pg/mL)	(pg/mL)	(pg/mL)	(pg/mL)	(pg/mL)	(pg/mL)
18–31 years (n = 41)	866.9	352.8	1752.18	56.05	137.69	58,537.07	638.57	229.37
18–31 years (n = 41)	(±1213.17)	(±446.45)	(±3198.11)	(±47.35)	(±144.44)	(±80,553.72)	(±1598.69)	(±223.48)
F (n = 24)	1008.46	369.69	1501.18	73.12	163.47	77,620.82	720.5	268.91
F (n = 24)	(±1488.28)	(±459.08)	(±2377.86)	(±53.41)	(±173.43)	(±99,722.87)	(±2023.03)	(±239.51)
M (n = 17)	667.06	328.97	2106.52	31.96	101.3	31,595.29	522.91	173.54
M (n = 17)	(±650.34)	(±440.81)	(±4148.96)	(±21.23)	(±80.72)	(±24,684.21)	(±694.31)	(±191.68)
32–44 years (n = 40)	798.12	246.3	1022.4	55.85	152.69	46,427.69	822.39	332.26
32–44 years (n = 40)	(±856.73)	(±213.42)	(±2056.48)	(±38.18)	(±289.03)	(±72,142.26)	(±985.87)	(±291.54)
F (n = 21)	1011.63	235.77	1301.05	65.26	200.19	65,836.55	929.3	339.12
F (n = 21)	(±1099.11)	(±225.48)	(±2268.67)	(±32.36)	(±391.17)	(±95,687.15)	(±833.59)	(±314.94)
M (n = 19)	562.13	257.95	714.42	45.45	100.18	24,975.79	704.22	324.68
M (n = 19)	(±368.87)	(±204.75)	(±1803.78)	(±42.16)	(±73.78)	(±13,335.03)	(±1142.67)	(±271.72)
45–62 years (n = 37)	419.26	261.82	2184.8	42.41	180.43	31,448.49	1842.39	237.24
45–62 years (n = 37)	(±429.22)	(±300.77)	(±7354.85)	(±23.44)	(±363.04)	(±21,726.21)	(±3308.21)	(±142.3)
F (n = 20)	247.93	175.63	894.66	45.77	136.21	36,615	1373.43	202.3
F (n = 20)	(±210.12)	(±274.52)	(±1392.53)	(±25.7)	(±125.51)	(±21,039.07)	(±2411.03)	(±111.15)
M (n = 17)	620.84	363.22	3702.61	38.47	232.47	25,370.24	2394.11	278.35
M (n = 17)	(±531.24)	(±306.32)	(±10,718.21)	(±20.52)	(±522.03)	(±21,536.14)	(±4138.03)	(±166.04)
63–88 years (n = 39)	478.46	272.72	1341.34	34.68	215.38	38,545.33	1170.57	207.61
63–88 years (n = 39)	(±579.94)	(±355.27)	(±3136.48)	(±21.45)	(±318.38)	(±48,499.9)	(±1112.74)	(±180.62)
F (n = 17)	310.61	248.91	1663.98	27.3	194.13	34,047.76	1239.38	153.99
F (n = 17)	(±379.82)	(±443.5)	(±4396.99)	(±19.87)	(±384.55)	(±48,595.44)	(±1158.81)	(±125.24)
M (n = 22)	608.16	291.11	1092.03	40.38	231.8	42,020.72	1117.4	249.04
M (n = 22)	(±676.95)	(±278.79)	(±1709.24)	(±21.3)	(±264.78)	(±49,276.32)	(±1100.28)	(±207.26)
Total (n = 157)	647.39	284.33	1566.15	47.48	170.88	44,101.84	1101.26	252.03
Total (n = 157)	(±850.46)	(±340.4)	(±4329.85)	(±35.53)	(±286.76)	(±61,226.66)	(±1983.23)	(±221.37)
F (n = 82)	679.1	263.02	1335.75	54.94	172.58	55,568.04	1040.8	246.82
F (n = 82)	(±1049.42)	(±366.63)	(±2688.73)	(±40.1)	(±282.62)	(±77,538.17)	(±1741.18)	(±226.78)
M (n = 75)	612.73	307.64	1818.05	39.32	169.03	31,565.46	1167.36	257.73
M (n = 75)	(±563)	(±309.98)	(±5611.22)	(±27.79)	(±293.12)	(±31,887.55)	(±2228.15)	(±216.69)

Pg/mL = picograms per milliliter. F = Female; M = Male; EGF = Epidermal Growth Factor; IL = Interleukin; RA = Receptor Antagonist; IP = Induced Protein; MCP = Monocyte Chemotactic Protein; VEGF = Vascular Endothelial Growth Factor. Mean (±Standard Deviation (SD)) is shown considering four equally size groups of age.

Table 3. Comparisons of fitted GAMLSS models.

	Normal		Box-Cox Cole and Green		Power Exponential		Box-Cox Power Exponential
	Model Specification	AIC	Model Specification	AIC	Model Specification	AIC	Model Specification	AIC
EGF	μ ~ age σ ~ 1	628.43	μ ~ age σ ~ 1 ν ~ age	629.3	μ ~ age σ ~ 1 ν ~ 1	629.9	μ ~ age σ ~ 1 ν ~ 1 τ ~ 1	627.38
Fractalkine	μ ~ sex σ ~ 1	576.18	μ ~ 1 σ ~ 1 ν ~ 1	579.37	μ ~ sex σ ~ 1 ν ~ 1	577.73	μ ~ 1 σ ~ 1 ν ~ 1 τ ~ 1	580.88
IL-1RA	μ ~ cs(age) σ ~ 1	690	μ ~ cs(age) σ ~ 1 ν ~ 1	687.66	μ ~ cs(age) σ ~ 1 ν ~ cs(age)	687.71	μ ~ cs(age) σ ~ 1 ν ~ 1 τ ~ 1	687.59
IL-7	μ ~ cs(age) + sex σ ~ 1	442.84	μ ~ age + sex σ ~ 1 ν ~ 1	445.65	μ ~ cs(age) + sex σ ~ age ν ~ 1	444.98	μ ~ age + sex σ ~ 1 ν ~ 1 τ ~ age + sex	442.20
IL-8	μ ~ cs(age) + sex σ ~ cs(age)	541.76	μ ~ cs(age) + sex σ ~ 1 ν ~ cs(age)	538.11	μ ~ cs(age) + sex σ ~ 1 ν ~ 1	531.89	μ ~ cs(age) + sex σ ~ 1 ν ~ 1 τ ~ 1	533.83
IP-10	μ ~ age + sex σ ~ age + sex	567.61	μ ~ sex σ ~ 1 ν ~ cs(age)	562.86	μ ~ age + sex σ ~ age + sex ν ~ 1	554.07	μ ~ sex σ ~ 1 ν ~ 1 τ ~ 1	556.26
MCP-1	μ ~ age σ ~ 1	672.06	μ ~ age σ ~ 1 ν ~ cs(age)	674.46	μ ~ age σ ~ 1 ν ~ cs(age)	672.32	μ ~ age σ ~ 1 ν ~ cs(age) τ ~ cs(age) + sex	673.6
VEGF	μ ~ cs(age) σ ~ 1	524.07	μ ~ cs(age) σ ~ 1 ν ~ 1	527.26	μ ~ cs(age) σ ~ 1 ν ~ age	517.01	μ ~ cs(age) σ ~ 1 ν ~ 1 τ ~ 1	521.21

Four different distribution families: Normal, Box-Cox Cole and Green, Power Exponential and Box-Cox power exponential were tested. Potential explicative variables are sex and age. The age influence on each distributional parameter (location: μ; scale: σ; skewness: ν; kurtosis: τ) can be modelled as a constant (~1), as a linear function (~age) or as a cubic spline (~cs(age)). Sex can be added as qualitative variable (~sex). The best model for a particular molecule, in bold, is the one that minimizes the Akaike Information Criterion (AIC). EGF = Epidermal Growth Factor; IL = Interleukin; RA = Receptor Antagonist; IP = Induced Protein; MCP = Monocyte Chemotactic Protein; VEGF = Vascular Endothelial Growth Factor.

Table 4. Validation and utility of reference intervals (RIs).

	Validation Sample				Utility Sample
	Valid n	Low	Normal	High	Valid n	Low	Normal	High
	Valid n	% (CI95%)	% (CI95%)	% (CI95%)	Valid n	% (CI95%)	% (CI95%)	% (CI95%)
EGF	40	0% (0%; 10.9%)	100% (89.1%; 100%)	0% (0%; 10.9%)	13	23.1% (6.2%; 54%)	76.9% (46%; 93.8%)	0% (0%; 28.3%)
Fractalkine	36	0% (0%; 12%)	94.4% (80%; 99%)	5.6% (1%; 20%)	11	0% (0%; 32.1%)	27.3% (7.3%; 60.7%)	72.7% (39.3%; 92.7%)
IL-1RA	39	0% (0%; 11.2%)	92.3% (78%; 98%)	7.7% (2%; 22%)	13	0% (0%; 28.3%)	7.7% (0.4%; 37.9%)	92.3% (62.1%; 99.6%)
IL-8	40	0% (0%; 10.9%)	95% (81.8%; 99.1%)	5% (0.9%; 18.2%)	13	0% (0%; 28.3%)	7.7% (0.4%; 37.9%)	92.3% (62.1%; 99.6%)
IP-10	35	2.9% (0.2%; 16.6%)	97.1% (83.4%; 99.9%)	0% (0%; 12.3%)	13	38.5% (15.1%; 67.7%)	61.5% (32.3%; 84.9%)	0% (0%; 28.3%)
MCP-1 *	8	0% (0%; 40.2%)	100% (59.8%; 100%)	0% (0%; 40.2%)	1	100%	0%	0%
VEGF	32	0% (0%; 13.3%)	93.8% (77.8%; 98.9%)	6.3% (1.1%; 22.2%)	12	0% (0%; 30.1%)	25% (6.7%; 57.2%)	75% (42.8%; 93.3%)

Two external samples are evaluated: (1) validation sample, consisting in the cytokine levels of 40 healthy subjects, and (2) utility sample that involved 13 severe cases of dry eye disease (DED). The percentage of values falling both within and outside the lower and upper reference limits is showed. CI = Confidence interval. * For MCP-1, there was only a valid subject in utility sample and he was classified as Low. IL-7 RI was not validated as its concentration information was missing in all subjects of both samples.

Publisher’s Note: MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

© 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https://creativecommons.org/licenses/by/4.0/).

Share and Cite

MDPI and ACS Style

Fernández, I.; Enríquez-de-Salamanca, A.; Portero, A.; García-Vázquez, C.; Calonge, M.; Herreras, J.M. Age- and Sex-Adjusted Reference Intervals in Tear Cytokine Levels in Healthy Subjects. Appl. Sci. 2021, 11, 8958. https://0-doi-org.brum.beds.ac.uk/10.3390/app11198958

AMA Style

Fernández I, Enríquez-de-Salamanca A, Portero A, García-Vázquez C, Calonge M, Herreras JM. Age- and Sex-Adjusted Reference Intervals in Tear Cytokine Levels in Healthy Subjects. Applied Sciences. 2021; 11(19):8958. https://0-doi-org.brum.beds.ac.uk/10.3390/app11198958

Chicago/Turabian Style

Fernández, Itziar, Amalia Enríquez-de-Salamanca, Alejandro Portero, Carmen García-Vázquez, Margarita Calonge, and José M. Herreras. 2021. "Age- and Sex-Adjusted Reference Intervals in Tear Cytokine Levels in Healthy Subjects" Applied Sciences 11, no. 19: 8958. https://0-doi-org.brum.beds.ac.uk/10.3390/app11198958

Note that from the first issue of 2016, this journal uses article numbers instead of page numbers. See further details here.

Article Menu

Age- and Sex-Adjusted Reference Intervals in Tear Cytokine Levels in Healthy Subjects

Abstract

Featured Application

Abstract

1. Introduction

2. Materials and Methods

2.1. Study Sample

2.2. Tear Collection and Tear Inflammatory Molecule Analysis

2.3. Data Analysis

2.3.1. Reference Intervals Determination

2.3.2. Model Diagnostics

2.3.3. Reference Intervals Utility

3. Results

4. Discussion

5. Conclusions

Author Contributions

Funding

Institutional Review Board Statement

Informed Consent Statement

Data Availability Statement

Conflicts of Interest

Appendix A. Residual Analysis

Appendix B. Centile Curves Estimation

Appendix C. Reference Intervals Classification

References

Share and Cite

Article Metrics

Article Access Statistics

Further Information

Guidelines

MDPI Initiatives

Follow MDPI