Loss of adaptive capacity in asthmatics revealed by biomarker fluctuation dynamics upon experimental rhinovirus challenge

Asthma is a dynamic disease, in which lung mechanical and inflammatory processes often interact in a complex, unpredictable manner. We hypothesize that this may be explained by respiratory disease-related systems instability and loss of adaptability to changing environmental conditions, resulting in highly fluctuating biomarkers and symptoms. Using time series of inflammatory (eosinophils, neutrophils, FeNO), clinical and lung function biomarkers (PEF, FVC, and FEV1), we estimated this loss of adaptive capacity (AC) during an experimental perturbation with a rhinovirus in 24 healthy and asthmatic volunteers. Loss of AC was estimated by comparing similarities between pre- and post-challenge time series. Unlike healthy participants, the asthmatic’s post-viral-challenge state resembled significantly more other rhinovirus-infected asthmatics than their own pre-viral-challenge state (hypergeometric-test: p=0.029). This reveals loss of AC, and supports the novel concept that not only single physiological mechanisms, but interacting dynamic disease properties are altered in asthma and contribute to a more vulnerable phenotype.


Introduction
The dynamics of (patho)physiology are particularly relevant in chronic diseases. In uncontrolled conditions, the patients are almost always clinically unstable, exhibiting a fluctuation pattern of symptoms, including unexpected loss of control or even sudden, potentially life-threatening exacerbations. Asthma, a chronic airway inflammation with variable airway obstruction is an archetypical example of such a dynamic disease. It is well known that the relationship between the strength of the environmental trigger, and the resulting degree of airway inflammation, airway obstruction, and symptoms is highly nonlinear, thus asthma is often insufficiently controlled with anti-inflammatory drugs (www.ginasthma.org) (1).
The temporal behavior of disease biomarkers in asthma is related to the inherently complex, interacting pathophysiology of the respiratory system, involving environmental, inflammatory and airway mechanical components (2). We and others have hypothesized that in health, the respiratory system exhibits a so-called 'homeokinetic' stable equilibrium state, enabling the respiratory system to easily adapt to an external perturbation (3). However, the homeokinesis of the respiratory system in asthma may be altered and adaptability to external stimuli (adaptive capacity) may be decreased, leading to a less stable system with more temporal fluctuations and exacerbations (3). Thus, asthma instability and fluctuation dynamics may not just relate to the degree of airway inflammation or the strength of the environmental trigger, but to the innate properties of the homeokinetic system itself.
Homeokinesis is the ability of a physiological system to return to a dynamic equilibrium after a perturbation (4-7); The remarkable complexity of such systems originates from non-linear interactions between their constitutive parts (8)(9)(10). The healthy homeokinetic respiratory system is characterized by a stable equilibrium and sufficient adaptability in the face of changing environmental conditions such as pollutants, pathogens and allergens (10). A concomitant phenomenon of homeokinesis is the fluctuation behavior of physiological processes and related physiological parameters and biomarkers around a stable medium value when measured longitudinally (11,12). Over the past decades, a considerable research effort has been invested into mathematically analyzing and characterizing physiologic time series (8,(13)(14)(15)(16). Moreover, it has been shown that time series properties, such as unpredictability, longrange correlations, fractality, and information content are associated with states of health and disease (10,17). Such changes in dynamic behavior are potentially related to loss of adaptive capacity of the complex respiratory system. However, quantifying such loss of the system's adaptive capacity, i.e., its ability to cope with external perturbations (10), is more difficult and dependent on the physiological context. For instance, researchers have directly linked the capacity of rats to adapt to environmental heat stresses to the ability of the animal's liver cells to rapidly express the heat shock protein HSP70 in high quantities (18). Other scientists have suggested "the capacity of a physiological system to bring itself autonomously back to the normal homeostatic range after a challenge" as a more workable definition of adaptive capacity (19). We found evidence for a loss of adaptive capacity by quantitatively comparing the similarities of means and fluctuations in the pre-and post-viral-challenge time series of biomarkers in an unassuming data driven manner, using hierarchical clustering (20) and statistically assessing the clusters found by means of enrichment analysis (21)(22)(23).
The aim of this study is to test the hypothesis, whether in asthmatics the adaptive capacity to a standardized environmental perturbation, such as an experimental viral challenge, is altered in comparison to healthy subjects. As a proof of concept, such changes in the homeokinetic system properties in disease would support evidence that not only single factors, such as inflammation, but also system properties contribute to disease dynamics and phenotype stability.
In a prospective, longitudinally designed study comprising healthy and asthmatic subjects, we measured time series of a set of standard lung functional and inflammatory/immune biomarkers two months prior to and one month following an experimental rhinovirus 16 (RV16) infection induced by controlled and deliberate inoculation of healthy and asthmatic volunteers. This choice was driven by the fact that rhinovirus (RV) infections in asthmatics have been found to be among the most prominent external triggers of acute worsening of asthma symptoms, and of asthma exacerbations and loss of control (24,25). Adaptive capacity to the rhinovirus infection was compared between asthmatic and healthy subjects for the abovementioned biomarkers.

Experimental rhinovirus challenge while monitoring cohort participants
In all cohort participants (12 healthy and 12 asthmatic volunteers), the biomarkers/parameters listed in Table 1 below were measured during two months before, and during one month immediately after deliberate experimental inoculation with rhinovirus, resulting in pre-and post-viral-challenge time series of each biomarker/parameter. Plots of the time series of each biomarker can be found in Supplementary Figures File 2. For the healthy and the asthmatics groups separately, summary statistics of the average before the viral challenge (average over 2 months) and after the viral challenge (average over 1 month) of each of these biomarkers/parameters can be found in the SI. Eosinophil and neutrophil cell density in nasal lavage fluid 1x weekly 3x weekly Table 1: Biomarkers/parameters measured in each cohort participant during two months before, and during one month immediately after deliberate experimental inoculation with rhinovirus. The corresponding sampling frequencies can be found in columns 2 and 3. See the Methods section below for more details on the study design, and on the measurement procedures and laboratory assays used. FEV1: forced expiratory volume in one second. FVC: forced vital capacity. PEF: peak expiratory flow. FeNO: fractional expired concentration of nitric oxide.

Hierarchical clustering of biomarker time series
In order to quantitatively establish the degree of similarity or "proximity" between two time series of a given biomarker, we used the Earth Mover's Distance (EMD), which regards each of the time series as a univariate empirical distribution of the biomarker at hand (see Methods and SI for more details). The pre-and post-challenge time series of individual biomarker time series from all participants (both healthy and asthmatics) were clustered using the EMD as the distance metric between the time series. The outcomes for the levels of exhaled nitric oxide (FeNO), and the percentage of eosinophils in nasal lavage fluid are presented here, whereas the results for the other biomarkers are presented in the SI.

i) Time series of exhaled Nitric Oxide (FeNO)
Findings are summarized in Table 2. The corresponding dendrogram is depicted in Figure 1. In brief, we found three clusters. Cluster 1 consists of four time series stemming from 2 asthmatics.
As can be read off of the dendrogram in Fig. 1 below, and of the distance matrix depicted in Panel C of Fig. 2 (see Methods section below), these two participants are prominently different from the rest (regarding their FeNO time series), and might be regarded as outliers. Cluster 2 contains more healthy participants than expected by chance. In other words, Cluster 2 is enriched in healthy participants. Conversely, due to the balanced design of the cohort (equal numbers of healthy and of asthmatic participants), Cluster 2 is also depleted of asthmatic participants, i.e., it contains fewer asthmatic participants than expected by chance. And finally, Cluster 3, which is enriched in asthmatic participants. The marked difference between the members of Cluster 1 and Cluster 3 regarding FeNO time series revealed by the hierarchical clustering evidences a higher heterogeneity among the asthmatics, as compared to the healthy participants. In Cluster 2, the tendency for infected participants to be clustered together with their corresponding uninfected counterpart is statistically significant (p-value = 0.038, see Table   2 below). This is not the case for Cluster 3. The difference in this regard between Cluster 2 (mainly healthy participants) and Cluster 3 (mainly asthmatic participants) is further underpinned by the fact that the distributions of cophenetic distances between the infected cluster members and their uninfected counterparts are statistically significantly different between these two clusters (p-value=0.033, see also Fig. 3 in the SI). More specifically, an average of 34.5 in Cluster 3, as opposed to 12.0 in Cluster 2.
In Table 3, the sub-clusters found within Clusters 2 and 3, respectively (marked with orange and blue rectangles in Figure 1), are analyzed in terms of enrichment in or depletion of pre-and post-challenge time series. This analysis provides evidence for a statistically significant separation of pre-and post-challenge time series within Cluster 3. Indeed, the union of subclusters 3.1 and 3.2 is enriched in pre-challenge time series (p-value=0.029, see Table 3 below), whereas subcluster 3.3 is enriched in post-challenge time series (p-value=0.029, see Table 3 below). Such a separation cannot be observed within Cluster 2.
A bootstrap based sensitivity analysis of these findings can be found in the SI.  Table 2: Composition, enrichment analysis, and grouping characteristics of the clusters found by comparison of each participant's pre-and post-challenge time series of FeNO. Enrichment is marked in bold letters, depletion in italics; the corresponding p-values were calculated using the hypergeometric test. The empirical p-values for the proportion of pre-and post-pairs were calculated using simulated permutations (see Methods section). A participant is fully represented in a given cluster if both their pre-and post-challenge time series of measurements are contained in the cluster. For example, the healthy participant "P08H" is fully represented in Cluster 2, as both their pre-and post-challenge time series of FeNO measurements are members of Cluster 2 (see Fig. 1 below). Partial representation corresponds to the scenario in which only one of the two time series (pre-and post-challenge) is a member of the cluster. For instance, the asthmatic participant " P07A" is only partially represented in Cluster 2, because their pre-challenge time series of FeNO measurements is part of Cluster 2, whereas their post-challenge time series of FeNO belongs to Cluster 3 (see Fig. 1 below). See also the Methods section for the definition of neighbors.

Figure 1:
Cluster dendrogram obtained via hierarchical clustering of the participants' pre-and post-challenge time series of FeNO. The distance between any two-time series was calculated using the EMD. Rectangles mark the clusters and sub-clusters identified. From top to bottom: Cluster 1, Cluster 2 (subdivided into Clusters 2.1, 2.2, and 2.3), and Cluster 3 (subdivided into Clusters 3.1 and 3.2, and 3.3). Patient IDs are indicated by Pxy, their health status using H/A, denoting Healthy or Asthmatic, and their RV infection status by Uninf/Inf, which stands for Uninfected/Infected. Cluster 1 consists of asthmatics which are prominently different from other asthmatic subjects in Cluster 3 and also from healthy subjects in Cluster 2. These might be regarded as outliers.

ii) Time series of percentage of eosinophils in nasal lavage fluid
Findings are summarized in Table 4. The corresponding dendrogram is depicted in Fig.1 in the SI. In brief, three clusters were identified. Cluster 1 consists of four time series stemming from 3 asthmatics; As can be read off of the dendrogram depicted in Fig. 1 in the SI, these time series are prominently different from all the other time series, and might be regarded as outliers.
Cluster 2 is enriched in healthy participants. And finally, Cluster 3, which is enriched in asthmatic participants. As seen in the analysis of FeNO, the marked difference between the members of Cluster 1 and Cluster 3 evidences a higher heterogeneity among the asthmatics, as compared to the healthy participants. However, Cluster 1 in the eosinophil analysis and Cluster 1 in the FeNO analysis only have one asthmatic patient in common. Again, in Cluster 2, the tendency for infected participants to be clustered together with their corresponding uninfected counterpart is statistically significant (p-value = 0.001, see Table 4 below). This is not the case for clusters 1 and 3. The difference in this regard between Cluster 2 (mainly healthy participants) and Cluster 3 (mainly asthmatic participants) is further substantiated by the fact that the distributions of cophenetic distances between the infected cluster members and their uninfected counterparts are statistically significantly different between these two clusters (p-  Table 4: Composition, enrichment analysis, and grouping characteristics of the clusters found by comparison of each participant's pre-and post-challenge time series of percentage of eosinophils in nasal lavage fluid. Enrichment is marked in bold letters, depletion in italics; the corresponding p-values were calculated using the hypergeometric test. The empirical p-values for the proportion of pre-and post-pairs were calculated using simulated permutations (see Methods section). A participant is fully represented in a given cluster if both their preand post-challenge time series of measurements are contained in the cluster. Partial representation corresponds to the scenario in which only one of the two time series (pre-and post-challenge) is a member of the cluster. See also the Methods section for the definition of neighbors.

Individual response to the viral challenge with respect to the biomarkers measured
In order to test the effectiveness of the virus challenge, we measured the individual patient's response with respect to each predefined biomarker. The results of blood antibody tests (RV16 seroconversion) along with RV Polymerase Chain Reaction (PCR) conducted on nasal lavage fluid taken from every participant after the inoculation indicated that 11 out of 12 healthy participants and 12 out of 12 asthmatics were effectively infected with the RV16 after inoculation ( Table 2 in the SI). According to these laboratory tests, one healthy participant did not become infected. However, this participant did develop cold symptoms within a few days after the virus inoculation, suggesting that the laboratory tests failed to detect the ongoing infection although the participant was positively infected. Consequently, this participant was included in the analyses.
After having established the efficacy of the inoculation with RV16, we then explored, for each of the biomarkers measured (listed in the first column of Table 5), for how many participants a statistically significant within-subject change upon infection can be observed ("responders", see Table 5). To this end, two criteria for "responders" were implemented. The first criterion,   Table 5: Proportions of responders within the groups of healthy and asthmatic participants, respectively. Two different criteria were used in order to establish a statistically significant response. According to the first criterion, a participant is considered a responder with respect to a given biomarker if the outcome of comparing the pre-challenge time series and the post-challenge time series of the same biomarker by means of the Kolmogorov-Smirnov test results in a p-value <= 0.05 (columns 2 and 3). According to the second criterion, a participant is considered a responder with respect to a given biomarker if the outcome of comparing, by means of a Mann-Whitney-U-test, the magnitude of relative changes observed during 10day time intervals prior to the challenge with the magnitude of relative changes that took place during 10-day time intervals that contained the day of the challenge results in a p-value <= 0.05 (columns 4 and 5). For calculating the proportion of responders within each group the p-values were corrected for multiple testing using the false discovery rate (FDR) method of Benjamini and Hochberg. FEV1: forced expiratory volume in one second. FVC: forced vital capacity. PEF: peak expiratory flow. FeNO: fractional expired concentration of nitric oxide. The lung function parameters FEV1 and FVC, and thereby their ratio FEV1/FVC, were normalized using the standardized reference equations recommended by Global Lung Function Initiative (GLI) Task Force for comparisons across different populations.

Discussion
In this proof of concept study, we provided experimental evidence for the loss of adaptive capacity in the human respiratory system due to asthma. To this end, we hypothesized that a loss of adaptive capacity could be experimentally demonstrated by means of a detection of a similarity diminution between the pre-and post-perturbation dynamics of the system. Using a data-driven clustering approach, we have shown that, in particular, FeNO and eosinophil time series were similar prior to and following the challenge in healthy subjects, suggesting stable homeokinetic behavior. In asthmatics, however, this similarity was predominantly reduced, suggesting a marked impact of the asthmatic condition on dynamic properties of the respiratory system, consistent with more unstable behavior and loss of adaptive capacity following the perturbation with viral infection. Adaptive capacity is a system characteristic of a homekinetic system. Our data on the respiratory system in asthma and health support the hypothesis that these homeokinetic system characteristics of the lung likely contribute to the stability and dynamic phenotype of asthmatic disease. This is particularly astonishing, since the average FeNO, lung function values, and nasal lavage fluid eosinophil counts are, in general, not different following a viral challenge.

The notion of adaptive capacity in a clinical and epidemiological context
Epidemiological asthma research has provided indirect clues about a loss of adaptive capacity in the human respiratory system due to asthma. For instance, respiratory comorbidities have been found to be remarkably more prevalent in asthma than in non-asthma patients (26).

Interpretation of our findings in a physiological context
Our study and findings for the first time demonstrate this at the system level using respiratory biomarker fluctuation dynamics. Loss of adaptive capacity cannot be determined by single-time point observations of biomarkers, it was, however, possible in a time series experiment involving asthma biomarkers during a standardized perturbation of the respiratory system. Nevertheless, the following considerations have to be discussed. We started off with the question: For which type of participant, healthy or asthmatic, and for which biomarker is the disruption introduced by the viral challenge strong enough to render infected individuals more similar among themselves than to their uninfected counterparts? Our results indicate that, with respect to cell density in nasal lavage fluid, there is a statistically significant separation of preand post-challenge time series. However, this effect does not seem to be cluster specific, even though we found a cluster of size 7 enriched in asthmatic participants within the subgroup of mainly post-challenge time series (see SI). A more distinctive scenario results from our cluster analysis of the pre-and post-challenge time series of the percentage of eosinophils in nasal lavage fluid: In Cluster 2, which is statistically significantly enriched in healthy participants, the tendency for infected participants to be clustered together with their corresponding uninfected counterpart is statistically significant and clearly higher than in Cluster 3, which is mainly composed of asthmatic participants. Finally, in the clustering of the pre-and postchallenge time series of FeNO we found Cluster 2, mainly composed of healthy participants, and Cluster 3, made of nearly 95% asthmatics. Furthermore, within Cluster 3 we found a statistically significant separation of pre-and post-challenge time series. No such separation was found within Cluster 2.

Effectiveness of the perturbation (viral challenge)
Our analyses were systematically carried out with the general goal of testing the hypothesis of a loss of adaptive capacity of the respiratory system in asthmatic cohort participants. However, towards this goal, it was pertinent to determine whether and in which individuals the viral inoculation was effective. The presence of cold symptoms and/or positive blood antibody tests along with RV PCR provide strong evidence of all participants being effectively challenged shortly after RV16 inoculation. Indeed, RV PCR was performed on day 3 post inoculation and in all but one participant the presence of RV was detected. This is in concordance with other studies showing that RV concentrations in the nasal wash could be detected 2-3 days after viral inoculation, with a signal peak on the 5 th /6 th day (39,40).

Variable and heterogeneous effect on lung function and inflammatory/immune biomarkers
We carried out a quantitative characterization of individual response to the viral perturbation.
This was done using two computational/statistical approaches. One approach aimed to capture the changes elicited by the viral challenge taking place over longer time periods (comparison of the per-and post-challenge time series, viewed as empirical distributions), whereas the other assessed relative short-term changes occurring at shorter time scales (comparison of the magnitude of relative changes observed during 10-day time intervals).
There is a clear macroscopic/functional manifestation of the kindling RV infection, as reflected at the level of distributional changes induced by the viral challenge on lung function parameters.
Indeed, with respect to this criterion, 50% or more statistically significant responders in each of the two groups (healthy and asthmatics) were found (see rows 1-4 in Table 5 above). Notably, significant differences found between pre-and post-challenge time series were, in general, not attributable to changes in the variance, as verified using Levene's test (results not shown).
Nevertheless, for most participants, the lung function parameters did not show shortterm/transient relative changes induced by the viral challenge that were statistically significantly different in magnitude from short-term changes observed during the pre-challenge phase (see columns 4 and 5, and rows 1-4 in Table 5 above, and Figures in Supplementary   Figures File 1). Taken together, these results suggest that the changes in lung function elicited by the viral challenge are, both for healthy and asthmatic participants, subtle, spread over comparatively longer time periods, and unlike a transient decline. This is in line with the results of previous studies (41)(42)(43)(44), which concluded that after RV challenge lung function in asthmatic subjects did change, but did not decline dramatically in comparison to the changes observed in healthy controls. In contrast, our analyses indicate that changes in the inflammatory or immunological biomarkers at the cellular or molecular level are short-term and transient in nature (see rows 5-8 in Table 5 above, and Figures in Supplementary Figures File 1).
Nevertheless, for these parameters fewer responders were found, when compared to the lung function parameters. However, our results also hint at a relatively short time scale of response of these inflammatory/immunological biomarkers. Thus, the sampling frequency used in this study may not entirely capture the rapidly changing magnitudes. The observed differences in the type of response between the lung function and the inflammatory/immunological biomarkers may be a manifestation of the interplay of different temporal and spatial scales.

Limitations of the study
Our findings need to be judged in light of the limitations of our study. One of the limitations is that we only included mild asthmatics for ethical reasons. However, more severe asthmatics are likely to be on corticosteroid treatment, which may potentially influence the inflammatory response, therefore introducing a confounding factor. Another potential shortcoming of this study is the relatively small sample size. However, this drawback is compensated for by the unprecedented high sampling frequency at which the participants were screened in our study.
In fact, the sample size was a compromise between the number of subjects and the number of assessment visits. Nevertheless, as discussed above, the sampling frequency for some of the inflammatory/immunological biomarkers may still not have been high enough in order to capture their temporal dynamics.

Conclusion and implications
We have found evidence supporting the notion that a chronic disease such as asthma may alter the properties of a homeokinetic physiologic system in a way that compromises its capacity to appropriately react to a possibly harmful environmental stimulus. This loss of adaptive capacity in the asthmatic lung may be understood as changes that render the system overly unstable (3).
As a proof of concept, such changes in homeokinetic system properties would support evidence that not only single factors, such as inflammation, but also system properties contribute to disease dynamics and phenotype stability. Future studies involving time series of biomarker measurements may help us understand this system instability in chronic asthma, even in the absence of severe airway inflammation or obstruction. Furthermore, future therapeutic approaches may want to focus on maintaining a stable homeokinetic equilibrium of the respiratory system, rather than just normalizing single physiological or inflammatory biomarkers.

Materials and Methods
This study was approved by the medical ethical committee from the Amsterdam University Medical Center and registered at the Netherlands Trial Register (NTR5426).

Participant cohort
Twelve non-smoking, atopic (as determined by positive skin prick test to common aero allergens), mild to moderate asthmatic subjects (based on ATS/ERS criteria), not using steroids were chosen for inclusion. Similarly, 12 non-smoking, non-atopic healthy subjects were also included in the study as controls. All participants provided written informed consent. The demographics of the study population are summarized in Table 7. All the participants were required to have their serum antibody titer of RV16 < 1:8 during screening. The age group for the study population was 18-30 years. Individuals with concomitant disease and pregnant women were excluded.
The basic inclusion criteria for the study populations followed standard recommendations and were as shown in Table 6 below.

Study Design
The project represents a prospective observational, follow-up study including patients with asthma and healthy controls with an experimental RV intervention.
The study participants were recruited after meticulous screening of volunteers (as mentioned in SI). The study was mainly divided into 2 phases. Phase 1 (stable phase) consisted of 2 months where these subjects were followed up and sampled diligently every alternate day (3 times a week for most of the measurements) with constant frequency at the hospital clinic. After that they were subjected to a standardized nasal dose of RV inoculation in the laboratory and followed up at the same frequency for one additional month (also called Phase 2 or unstable phase). In total the study consisted of 3 months of sampling period with a minimum of 180 measurements of lung function, 33 FeNO data points, and 20 cytokine and cell count measurements per subject.
The schematic work flow of the phases mentioned, is provided in the Fig 10 of the SI.

Lung function assessment
Spirometry was performed only once on the screening visit at the clinic to include participants based on inclusion criteria using a daily calibrated spirometer according to European Respiratory Society (ERS) recommendations (45).
Home monitoring of morning and evening lung function was done by hand held devices (Micro Diary, CareFusion, yielding the FEV1, FVC, FEV1/FVC and PEF values analyzed in this study.
Moreover, the Asthma Control Questionnaire was administered.

Exhaled Nitric Oxide (FeNO)
Measurement of fractionated exhaled nitric oxide (FENO) was performed using the NIOX MINO® (Aerocrine AB, Sweden). Single measurements per person were recorded at the clinic, thrice weekly, according to recommendations by the ATS (46).

Nasal lavage
Nasal lavage was collected from the study participants once weekly before RV challenge and was up scaled to thrice weekly after the challenge at the clinic as previously described (42) [Refer to SI for details].

Rhinovirus challenge
The study participants were exposed to rhinovirus 16 (RV16) using a standardized and validated challenge approach, based on previous studies by ourselves and other groups (47)(48)(49)(50)(51)(52) show that a low dose is sufficient to induce mild cold-symptoms (52,53). Furthermore, this low-dose inoculum previously resulted in a slight decrease of FEV1 (loss of asthma control) in asthmatic patients between day 4 and 6 after RV16 exposure, whereas no decrease has been observed in healthy controls (50).
Refer to SI for further details.

Statistical and computational Analysis
Statistical tests resulting in a p-value less or equal to 0.05 were regarded as significant. Our clustering approach makes use of the entire time series (distributions) of values measured before and after the challenge, respectively, and does not amalgamate the information into a single magnitude (e.g., the mean value). This method unveils subtle differences and similarities between the participants' measurements that are less likely to be captured by conventional methods based on averages.  The start position of the gliding window is expressed relative to the day of the viral challenge, which is marked as day 0. When the position of the gliding window was such that the day of the viral challenge was contained within the gliding window, the corresponding value of the relative change is marked in red. In order to assess the statistical significance of the short-term relative changes possibly elicited by the viral challenge, the relative change values located to the left of those marked in red were compared to the values marked in red by means of a Mann-Whitney-U-test. Visual inspection of the time series in Panel B correctly suggests that the outcome of this test is not significant. The reason being that the relative changes within time intervals of 10 days observed prior to the viral challenge are comparable to changes observed within intervals of the same length containing the day of the viral challenge. In Panel C, depicting data from a different healthy participant, the situation is clearly different, as verified by a significant outcome of the corresponding Mann-Whitney-U-test. In such cases, the participant is called a "responder" with respect to the "relative change within 10 days criterion".

Characterization of the dendrogram clusters
In order to evaluate the discriminatory power of a given biomarker, the clusters found in the clustering dendrogram were tested for enrichment in or depletion of healthy or asthmatic participants, and/or for enrichment in or depletion of pre-or post-challenge distributions.
The relative location of leaves in the clustering dendrogram was quantitatively evaluated using the cophenetic distance (58). The cophenetic distance between two leaves of a dendrogram is defined as the height of the dendrogram at which the two largest branches that individually contain the two leaves merge into a single branch.
For every cohort participant and any given biomarker there is a pre-challenge and a postchallenge time series, which we call the participant's pre-and post-pair. If the disruption caused by the viral challenge is not strong enough, the pre-and post-challenge distributions of a given participant will tend to cluster together. Therefore, a cluster in which pre-and post-pairs are closely located in terms of the cophenetic distance within the dendrogram, represents a subgroup of participants for which the viral challenge caused a relatively weaker disruption, at least with respect to the biomarker under scrutiny.
Two dendrogram leaves are called neighbors if their mutual cophenetic distance is equal to the minimum of all cophenetic distances from one of the leaves to all the other leaves in the dendrogram. If this condition is fulfilled for both leaves simultaneously, then the two leaves form a two-element cluster in the dendrogram. If the condition is only fulfilled for one of the leaves, the two are still considered neighbors, even if this is not always visually obvious from inspecting the dendrogram (see Figure 11 in the SI).
Under the null-hypothesis that the branching in the dendrogram is the result of a purely random process, the number of neighboring pre-and post-pairs to be expected just by chance within a given cluster can be estimated by simply permuting the labels of the leaves in the dendrogram and counting the number of neighboring pre-and post-pairs. This permutation test is used for calculating the empirical p-values displayed in Tables 2 and 4

Effectiveness of the viral inoculation:
Each participant in the study was administered the same dose of the virus (100 TCD 50) through the nose and every subject was tested for being positive for the virus after inoculation. False positive results due to previous exposure to the virus was ruled out by strict inclusion criteria of not having the titer of antibodies against RV16 > 1:8 in serum, measured at screening and prior to inoculation.
Positivity to viral inoculation was confirmed by either one of these three criteria

1) Positive test for antibodies against RV at the terminal visits of the participant 2) Positive RV PCR test from approximately 3 rd day (second visit) post RV challenge 3) Symptoms of RV induced cold
The effect of the virus in every study participant is summarized below in Table 1 (Column 1 indicates the seroconversion of the antibodies against the virus at the end of the study, Column 2 reflects the response to the virus by the PCR product on the day 3 after challenge and finally Column 3 shows the symptoms developed in the volunteers after the viral inoculation). Positive response in either of the 3 categories was considered evidence for a successful viral challenge.

Percentage of eosinophils in nasal lavage fluid
The corresponding dendrogram is depicted in Fig. 1 below.   Figure 3: The boxplot to the left represents the distribution of cophenetic distances between time series corresponding to infected healthy participants and their uninfected counterparts. Only time series belonging to Cluster 2 in the clustering dendrogram obtained using FeNo data (see Fig. 1 in the Main Manuscript) are contemplated here. The boxplot to the right represents the distribution of cophenetic distances between time series corresponding to infected asthmatic participants and their uninfected counterparts. Only time series belonging to Cluster 3 in the clustering dendrogram obtained using FeNo data (see Fig. 1 in the Main Manuscript) are contemplated here. The two distributions are statistically significantly different (p-value=0.033, one-tailed Mann-Whitney-U-test, the cophenetic distances in Cluster 3 being, on average, higher than the ones in Cluster 2.

Cell density in nasal lavage fluid:
We

Percentage of neutrophils in nasal lavage fluid
We found one cluster of size 10 enriched in time series from healthy participants (p= 0.036).
Moreover, no clusters were found showing a significant enrichment in or depletion of prechallenge or post-challenge time series.

Sensitivity analysis of main clustering results:
We investigated the sensitivity of the clustering of FeNO time series to changes in the data via nonparametric bootstrapping. However, given that the post-viral challenge time series are very short (eleven data points or fewer), resampling would be strongly affected by small sample size effects [1]. Thus, we only applied bootstrapping to the pre-challenge time series. Moreover, we resorted to soft bootstrapping [2] (see Methods below for more details) in order to increase the likelihood that least frequent values in the time series would be chosen during the resampling procedure. The results from 1000 soft bootstrapping iterations were as follows: In 100% of the resulting soft bootstrap dendrograms, Cluster 1 ("outliers" cluster) was found.

Participant cohort
24 study subjects were recruited among which 12 were healthy volunteers and the other 12 steroids naïve (or stopped using steroids 6 weeks prior to the study) mild to moderately persistent asthmatics.
The inclusion and exclusion criteria for participants to enter the study were as follows:

Inclusion criteria
Asthma patients were selected using the following inclusion criteria: • Age 18-50 years • History of episodic chest tightness and wheezing • Intermittent or mild to moderate persistent asthma according to the criteria by the Global Initiative for Asthma (Global Initiative of Asthma. www.ginasthma.org) • Non-smoking or stopped smoking more than 12 months ago and 5 pack years or less • Clinically stable, no exacerbations within last six weeks prior to study • Steroid-naïve or those participants who are currently not on corticosteroids and have not taken any corticosteroids by any dosing-routes within 6 weeks prior to the study or only using on-demand reliever therapy • Baseline pre-bronchodilator FEV1 ≥ 70% of predicted [3] • Airway hyperresponsiveness, indicated by a positive acetyl-ß-methylcholine bromide (MeBr)challenge with PC20£ 9.8 mg/ml [4] • Positive skin prick test (SPT) to one or more of the 12 common aeroallergen extracts, defined as a wheal with an average diameter of ³ 3 mm • No other clinically significant abnormality on history and clinical examination • Able to give written and dated informed consent prior to any study-specific procedures Healthy subjects were selected using the following inclusion criteria • Age 18-50 years • Non-smoking or stopped smoking more than 12 months ago and 5 pack years or less.
Steroid-naïve, non-atopic participants who are currently not on any maintenance (subjects using oral contraceptives can be accepted) • Baseline FEV1 ≥ 80% of predicted • Negative acetyl-ß-methylcholine bromide (MeBr) challenge or PC20 ≥ 19.6 mg/ml • Negative skin prick test (SPT) to one or more of the 12 common aeroallergen extracts • Negative history of pulmonary and any other relevant disease • Able to give written and dated informed consent prior to any study-specific procedure.

Exclusion criteria
Potential subjects who meet any of the following criteria were excluded from participation in the study: • Women who are pregnant, lactating or have a positive urine pregnancy test at baseline visit • Participation in any clinical investigational drug treatment protocol within the preceding 5 half-lives (or 12 weeks, if the half-life is unknown) before the screening visit of this study.
• Concomitant disease or condition which could interfere with the conduct of the study, or for which the treatment might interfere with the conduct of the study, or which would, in the opinion of the investigator, pose an unacceptable risk to the patient Furthermore, the following additional exclusion criteria were used in phase 2 of the study: • RV16 titer > 1:8 in serum, measured at screening (visit 1) and also at the visit prior to the challenge • History of clinically significant hypotensive episodes or symptoms of fainting, dizziness, or light-headedness • History of an asthma exacerbation within the last 6 weeks prior to the study • Has had any acute illness, including a common cold, within 4 weeks prior to visit 1 • Close contact with young children or with any immunosuppressed patients • Has donated blood or has had a blood loss of more than 450 mL within 60 days prior to screening visit or plans to donate blood during the study • Positive for any virus in nasal lavage at visit immediately prior to rhinovirus challenge.

Sample size calculation
The study is explorative in nature and is based on some first estimates of fluctuating inflammatory biomarkers. A definite sample size for all biomarkers measured was not possible to be accurately determined.
A sample size of 12/12, based on previous studies with lesser data points, was calculated to identify biomarkers that have a strong, clinically relevant, sensitivity to the changes of disease stability over time [5,6]. The sample size will probably have missed effects on biomarkers with weak effects, nevertheless, it has for all biomarkers screened provided evidence on long time variability, correlation and cross-correlation between biomarker properties, and susceptibility to strong viral external challenges. These three dimensions are a first estimate to understand the degree of variability and complexity of any disease process.

Study Design
The project represents a prospective observational, follow up study including patients with asthma and healthy controls with an experimental rhinovirus intervention in between.

Screening Visit
First, the subjects were asked to sign an informed consent form after which they were examined for inclusion and exclusion criteria along with medical history and scoring of adverse events and concomitant medication. Spirometry, Methacholine challenge and Skin Prick Test along with physical examinations were performed meticulously.

Run-in Phase
During this phase adverse events and medication were recorded along with physical examination. Exhaled nitric oxide (FeNO), nasal lavage, measurements were performed to acquaint the study subjects with the study procedures. Starting from the run-in, once-daily home-monitoring by Asthma Control Questionnaires, symptom scores and twice daily maneuvers by the pocket spirometer were executed.

Baseline visit
The medication history (also adverse events and concomitant medication) was carefully recorded and thereby routine FeNO, and nasal lavage were measured. Pregnancy test was also performed.

Study period
• Phase 1: The study participants were asked to visit the hospital clinic thrice weekly for the aforementioned measurements on Monday, Wednesday and Fridays, during the first 60 days.
• Phase 2: The participants were again followed up for the same measurements for the next 30 days with the same sampling frequency after being inoculated nasally with the common cold virus Rhinovirus 16, thereby mimicking the trigger for loss of asthma control or a mild exacerbation.
All participants in the study were screened for the presence of respiratory viruses just before the rhinovirus challenge, using Polymerase Chain Reaction confirmation on nasal lavage samples. This was performed in order to rule out a concomitant infection of the respiratory airways with viruses other than the one used in our study.

End of study visit
All measurements were repeated at the end of each phase.
The schematic representation of the phases of the study is mentioned below in Fig. 10

Asthma Control Questionnaire
The Asthma Control Questionnaire (ACQ) was used in this study to assess the symptoms before and after Rhinovirus challenge along with lung function coupled to a hand-held spirometer device. Though widely used and well validated, it covers a 7-day time span and may not be optimal for a challenge protocol wherein symptoms change daily. Hence ACQ was used on a daily basis to record the daily changes in symptoms. Subjects complete the diary on rising in the morning and retiring at bedtime.

Skin prick test
Allergy skin prick tests were performed based on the position paper by the European Academy of Allergology and Clinical Immunology (EAACI) [7]. The test was done by placing a drop of each of the 12 solutions containing common aeroallergens on the skin, followed by needle pricks by small SPT-lancets. Histamine was used as a positive control and saline solution as a negative control. The test result was considered positive if the skin develops a red, itchy area, called a 'wheal', with an average diameter of at least 3 mm, 15 minutes after the prick. The outline of all wheals was marked with a water-soluble marker and transferred to a test result from using adhesive tape, so as to facilitate measuring of the diameter of the wheal. All test result forms were archived. If the patient has been using anti-allergic medications, they were instructed to stop 3 days prior to the test in the clinic.

Spirometry
Spirometry was performed using a daily calibrated spirometer according to European Respiratory Society (ERS) recommendations. Spirometry was done only once on the screening visit at the clinic. There was no use of bronchodilators (pre/post measurements) in the study.
The subjects were asked not to use reliever medications at least 6 hours before the test starts in the clinic. An experienced lung function analyst from the AMC performed the lung function tests throughout the study to reach optimal performance and to enhance reproducibility. Waking and bedtime FEV1 separately are of interest due to the exaggerated diurnal airflow variation seen in asthmatics. Subjects were encouraged to measure FEV1 at consistent times (upon waking between approximately 6am-9am and prior to bedtime between approximately 8pm-11pm, respectively -always prior to any bronchodilator therapy).

PEF: Peak Expiratory Flow
All the home recordings were monitored regularly (every 2 weeks) for uniformity and consistency of data. Stringent quality control of home lung function measurements was done at each visit.

Exhaled nitric oxide
Measurements of fractionated exhaled nitric oxide (FENO) were performed (one measurement per subject) using the NIOX MINO ® (Aerocrine AB, Sweden). This analyser has an accuracy of ± 5 ppb or max 10%. The precision is < 3 ppb of a measured value < 30 ppb and ± 10 % of measured value > 30 ppb. Every individual subject was measured using a single measurement as standardised in the hospital.
Subjects were asked to perform a slow vital capacity manoeuvre with a standardised expiratory flow of 50 ml/sec for as long as possible. Positive mouth pressure was applied to close the velum and prevent contamination with NO from nose and paranasal sinuses. Expired gas was sampled continuously at the mouth piece and mean FENO values at the end expiratory plateau were calculated.

Methacholine challenge
Methacholine is a cholinergic synthetic analogue of acetylcholine and acts directly on the airway smooth muscle resulting in bronchoconstriction. Measurement of airway responsiveness can be done using incremental inhaled doses of methacholine as bronchoprovocation test. This challenge test was performed using MeBr (acetyl-β-methylcholine bromide) according to the standardized tidal volume method that is operative in routine clinical diagnostics in the hospital [4].

Nasal lavage
Nasal lavage was collected from the study participants once weekly before rhinovirus challenge and was up scaled to thrice weekly after the challenge at the clinic.
Eight ml of saline solution was introduced into the nasal cavity by a catheter and maintained for 10 min before recovery, followed by filtration and removal of mucous and cells by centrifugation. Standardized washings collected from the nose was used for biomarker analyses.

Cell and slide preparation:
A cell suspension of 200,000 cells/ml were prepared in PBS. 2 sets of Cytospin apparatus (Thermo Shandon Cytospin 4) were assembled and the cytospin filter was pre wet by spinning at 550 rpm for 1 minute with 50 µl PBS. 100 µl of the cell suspension was added to each slide and centrifuged at 450 rpm for 2 min. The quality of slides and the density of the cells were evaluated using a phase contrast microscopy. If the cell density was too high or low, more slides were prepared using an adjusted volume of the cell suspension. A target of 4 good slides with no overlapping or clumping of cells was set to improve the overall quality and consistency.
If the total cell count (for both slides together i.e. less than 7500 cells/slide) was less than or equal to 15000 cells, then the sample was centrifuged and re-suspended in 100 µl and was loaded onto a single block.
The prepared slides were stained as soon as possible with Diff Quick 1 stain (Medion Diagnostics), dried and mounted on Depex.

Nasal cell cytology:
A differential cell count on a maximum of 400 inflammatory cells were performed. Numbers of eosinophils, neutrophils, macrophages/monocytes, lymphocytes, columnar and squamous epithelial cells were recorded. A report was prepared with the overall counts, the quality of cytospin and additional aspects such as the presence of mucous, cell debris, inclusions, eosinophil granules etc.

Blood venapunction:
15 ml of venous blood was collected at each scheduled visit (once weekly) to determine whether circulating antibodies against RV16 are present, for hematology, and other immunological biomarker assays. Collection of blood was done in standard EDTA and non-EDTA tubes for specific purposes. 10 ml of it was immediately be followed by centrifugation to obtain plasma [2000g for 10 minutes (min) at room temperature, RT] aiding in removal of RBCs and WBCs.
The supernatants were aliquoted for biomarker hunt.

Rhinovirus challenge
In this study we exposed the volunteers to a mild dose (100 TCID 50) of RV16 using a validated approach which has been previously shown to be sufficient in inducing mild cold symptoms and decrease in lung function.
We used the GMP RV16 stock that has been tested previously by the medical ethical commission at the hospital in AUMC, and also by the U-BIOPRED showing efficacy in terms of cold symptoms and viral replication at 100 TCID50, which is part of the accompanying IMPD. This GMP RV16 stock has been prepared from a seed virus in extensively characterized human volunteers as described in METC 2010_310 and was expanded and aliquoted by Charles River Laboratories (USA) under GMP conditions. This preparation was considered safe for in vivo testing in human volunteers during a scientific advice meeting at BFarM (Bonn, Germany, April 30, 2013). This testing was carried out as per the approved protocol from the Amsterdam UMC medical ethical committee.
Therefore, this part of our methodology is based on the distributional properties of each participant's measurements and neglects the time dimension.

Clustering approach
In our approach, for each biomarker, the measurements collected before and after the viral challenge were used to construct the individual empirical distribution of measurements for each study participant, before and after the challenge, respectively.
We performed agglomerative hierarchical clustering [22] of the aforementioned empirical distributions of biomarkers. Within the hierarchical clustering algorithm, the distances between the distributions were calculated using the Earth Mover's Distance [23], and the agglomeration procedure was done according to Ward's minimum variance method [22]. Intuitively speaking, the Earth Mover's Distance contemplates the pair of distributions to be compared as piles of sand and measures the effort that it would take to shovel one distribution into the shape and position of the other (see below for a more detailed description of this method). For each participant, the day of the challenge (i.e., the day of the inoculation) was labeled as "day 0". All days between recruitment and the day of the challenge were labeled with negative integer numbers. All days between the day of the challenge and the participant's final visit were labeled with positive integer numbers. Given that measurements were not conducted every day on a given participant (see study design above), for each biomarker considered in this study and for each participant separately, interpolation was used in order to have, for every biomarker one value for every day. Except for the total cell count in nasal lavage fluid, and the percentages of eosinophils and neutrophils in nasal lavage fluid, which were linearly interpolated, all biomarker time series were interpolated using cubic splines with natural boundary conditions (see, e.g., [24]). These interpolated values were only used for the assessment of the shortterm/transient response.

Calculation of short-term/transient changes
The time interval of 10 days was decided based on prior knowledge where differences in signals from biomarkers (lung function etc.) usually take not less than a week to subside As a sensitivity analysis, the same calculations were repeated using a 7-day time interval (data not shown). The outcomes were very similar to the ones obtained using the 10-day interval.
Consequently, a 10-day time period was considered optimal to calculate the short term/transient changes.

Characterization of the dendrogram clusters
The clusters obtained using the clustering dendrogram were tested for enrichment in or depletion of healthy or asthmatic participants, and/or for enrichment in or depletion of pre-or post-challenge distributions.
The relative location of leaves in the clustering dendrogram was quantitatively evaluated using the cophenetic distance [25]. The cophenetic distance between two leaves of a dendrogram is defined as the height of the dendrogram at which the two largest branches that individually contain the two leaves merge into a single branch.
For every cohort participant and any given biomarker there is a pre-challenge and a post-challenge time series, which we call the participant's pre-and post-pair. If the disruption caused by the viral challenge is not strong enough, the pre-and post-challenge distributions of a given participant will tend to cluster together. Therefore, a cluster in which pre-and postpairs are closely located in terms of the cophenetic distance within the dendrogram, represents a subgroup of participants for which the viral challenge caused a relatively weaker disruption, at least with respect to the biomarker under scrutiny.
Two dendrogram leaves are called neighbors if their mutual cophenetic distance is equal to the minimum of all cophenetic distances from one of the leaves to all the other leaves in the dendrogram. If this condition is fulfilled for both leaves simultaneously, then the two leaves form a two-element cluster in the dendrogram. If the condition is only fulfilled for one of the leaves, the two are still considered neighbors, even if this is not always visually obvious from inspecting the dendrogram (see Fig. 11 below).
Under the null-hypothesis that the branching in the dendrogram is the result of a purely random process, the number of neighboring pre-and post-pairs to be expected just by chance within a given cluster can be estimated by simply permuting the labels of the leaves in the dendrogram and counting the number of neighboring pre-and post-pairs. This permutation test is used for calculating the empirical p-values displayed in Tables 1 and 3 in the Main Manuscript. Rectangles mark the clusters identified. Panel B displays a more detailed view of the second cluster. According to the definition of neighboring leaves provided in the text above, the leaves P03H_Uninf and P03H_Inf are neighbors in this dendrogram (highlighted in red). The reason for this is that the cophenetic distance from leaf P03H_Uninf to leaf P03H_Inf is equal to the minimum of all distances from leaf P03H_Uninf to all other leaves in the dendrogram. This is not the case for the leaf P03H_Inf. However, the fact that this condition holds for at least one of the two leaves renders them neighboring. Leaves P08H_Inf and P08H_Uninf, and P11H_Inf and P11H_Uninf, are neighbors, respectively (the latter pair is not marked in color). P08A_Inf and P08A_Uninf are also neighbors; In this case, the minimum condition is fulfilled by both leaves for of this leaf-pair. This is why the two leaves form a two-element cluster in the dendrogram. Panel C displays a more detailed view of the third cluster. Analogous information about neighboring leaves in Cluster 3 is highlighted in blue. However, as opposed to Cluster 2 (depicted in Panel B), the amount of neighboring leaves in Cluster 3 is not statistically significant (permutation test, see Table 2 in the Main Manuscript).

Soft Bootstrapping
As elucidated in [2], when resampling with replacement from a given time series in order to generate a bootstrap replicate of the same time series, the relative frequencies of the values in the original time series were adjusted in the sense of soft bootstrapping as follows: 1) The relative frequencies of the values in the original time series were sorted in ascending order: " ≤ $ ≤ % ≤ ⋯ ≤ ' , where m is the number of different values in the time series.

The Earth Mover's Distance
The Earth Mover's Distance [23] is a method for quantifying the dissimilarity between two probability distributions. Intuitively speaking, the EMD contemplates the pair of distributions to be compared as piles of sand and measures the minimal effort that it would take to shovel one distribution into the shape and position of the other.
In practice, two distributions will be given by two representative samples, which can be written as lists of pairs {(v1 , w1),…, (vn , wn)} and {(c1 , f1),…, (cn , fm)}. Each pair (vi , wi) corresponds to a value vi and its relative frequency wi in the sample. If we translate the above described intuitive approach into numbers, the problem becomes a well-known transportation problem [26]: Suppose that n suppliers are located at positions v1 ,…, vn, respectively, and each one has a given amount of goods wi. Furthermore, they are required to supply m consumers, located at positions c1 ,…, cm, respectively, whereas each one has a given specific demand fi. For each supplier-consumer pair (vi , wi) and (cj , fj), the cost of transporting a single unit of goods is determined by the distance d(vi , cj) between their locations. The transportation problem is then to find a least expensive pattern of flow of goods from suppliers to consumers that would satisfy the consumers' demand. Once the optimal pattern of the goods' flow has been found, the total cost is the corresponding EMD.
Mathematically, this transportation problem can be formalized as a linear programming problem, for which efficient solution algorithms were developed in the late 1940s (see, e.g. [27]).

Utility of biomarker time series analysis
Do longitudinal measurements provide deeper insights into complex disease physiology as compared to single measurements or average values? In order to answer this question, we tried to reproduce the results obtained using each participant's entire collections of pre-and postchallenge biomarker measurements after collapsing them to the corresponding pre-and postchallenge individual average.
For example, in this cohort, FeNO time series have the ability to discriminate between healthy and asthmatics, and, within the group of asthmatics, between infected and uninfected (see Fig.   1 in the Main Manuscript). In order to investigate whether the average value would have a similar discriminative power, we calculated, for each participant, the average of their pre-and of their post-challenge series and used the absolute value of the difference between averages as a distance measure for clustering. The resulting dendrogram is depicted in Fig. 12 below. While the discriminative power between healthy and asthmatics is still given, the ability to distinguish between infected and uninfected within the group of asthmatics gets lost.