Loss of adaptive capacity in asthmatic patients revealed by biomarker fluctuation dynamics after rhinovirus challenge

Asthma is a dynamic disease, in which lung mechanical and inflammatory processes interact in a complex manner, often resulting in exaggerated physiological, in particular, inflammatory responses to exogenous triggers. We hypothesize that this may be explained by respiratory disease-related systems instability and loss of adaptability to changing environmental conditions, manifested in highly fluctuating biomarkers and symptoms. Using time series of inflammatory (eosinophils, neutrophils, FeNO), clinical and lung function biomarkers (PEF, FVC,FEV1), we estimated this loss of adaptive capacity (AC) during an experimental rhinovirus infection in 24 healthy and asthmatic human 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 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 concept that in asthma, biological processes underlying inflammatory and physiological responses are unstable, contributing to loss of control.


Introduction
The quantitative study of physiologic systems, such as the respiratory system, has revealed their ability to maintain a highly organized internal environment that is fluctuating within certain limits, despite being constantly exposed to a variable external environment (Que et al., 2001;Tirone and Brunicardi, 2001;Palmer and Clegg, 2016). The term homeokinesis has been coined to describe this ability, substituting the concept of homeostasis (Goldstein and Kopin, 2017) to emphasize that fluctuations in the internal environment are normal (Macklem, 2008;Macklem and Seely, 2010;Yates, 1982;Eke et al., 2002;Glass, 2001). Homeokinesis is at the core of the observed adaptability, that is, adaptive capacity, of physiologic systems in response to changing environmental conditions (Goldberger et al., 2002). It is part of the remarkable complexity characteristic of such systems, which is believed to originate from non-linear interactions and feedback-loops between their constitutive parts (Goldberger et al., 2002;Garfinkel, 1983;Goldberger, 1996;Que, 1998).
Over the past decades a considerable research effort has been invested into mathematically analyzing the fluctuation behavior of physiologic time series with the aim of characterizing the normal homeokinetic variability of physiologic systems (Garfinkel, 1983;Glass and Mackey, 1988;Goldberger et al., 2000;Glass and Kaplan, 1993;Amigó and Small, 2017). One of the most counterintuitive findings has been that both excessive and too little variation are indicative of pathological modifications and aging (Que et al., 2001;Goldberger et al., 2002 and references therein). Within this paradigm, a chronic disease, such as asthma, may be understood as changes in the system that render it either too rigid or overly unstable (Frey et al., 2011). Consequently, such disease or aging related changes of the system are accompanied by a loss of adaptive capacity (Goldberger et al., 2002). However, whether this can be detected in asthmatics using longitudinal measurements, that is time series, of asthma-related biomarkers has never been investigated.
The aim of this study is to test whether in asthmatics the adaptive capacity to a standardized environmental perturbation, such as an experimental viral challenge, is altered in comparison to healthy subjects.
In this 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, asthma exacerbations, and of loss of control (Johnston et al., 1995;Nicholson et al., 1993).
Quantifying a loss of adaptive capacity, that is, an impairment in the ability to cope with external perturbations, in a physiologic system will depend on how adaptive capacity is defined. This, in turn, depends 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 (Hall et al., 2000). 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 (Oresˇic and Vidal-Puig, 2014). While still very general, the latter definition seems more suitable in a (patho-) physiological context in which the specific molecular mechanisms behind the adaptation processes cannot be easily laid out. Indeed, in geriatric medicine this definition has been widely used, and the term homeostenosis was introduced to describe aging related loss of adaptive capacity (Taffett, 2003;Fossion et al., 2018;Yashin et al., 2007;Fried et al., 2009).
Based on these ideas, we define adaptive capacity as the ability of a physiological system to autonomously return to the normal homeokinetic range after an external challenge. We implemented this definition quantitatively in our study by comparing the participants' pre-and post-viralchallenge time series of measurements of the aforementioned biomarkers. We hypothesized that, for a given biomarker, the post-viral-challenge time series of a study participant with an unimpaired adaptive capacity would resemble the same subject's pre-viral-challenge time series. Conversely, the post-viral-challenge time series of a study participant with an impaired adaptive capacity would be relatively distinct from the same subject's pre-viral-challenge time series. In order to test this hypothesis, hierarchical clustering was used to group time series according to their relative similarity in an unassuming, data-driven manner.
With this approach we found experimental evidence for the loss of adaptive capacity of the human respiratory system due to asthma.

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 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 Appendix.

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 Materials and methods and Appendix for more details). The pre-and post-challenge time series (also referred to as uninfected participant and infected participant, respectively) 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 Appendix.

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 two asthmatics. 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 Materials and 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.  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 Materials and 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 Figure 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 Figure 1 below). See also the Materials and methods section for the definition of neighbors. As can be read off of the dendrogram in Figure 1 below, and of the distance matrix depicted in Panel C of Figure 2 (see Materials and 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, that is it contains fewer asthmatic participants than expected by chance. And finally, Cluster 3, which is enriched in asthmatic participants. While all, but one, of the time series from healthy participants are grouped together in Cluster 2, the vast majority of time series from asthmatic participants are split into two different Clusters, namely Clusters 1 and 3. This suggest a higher heterogeneity among the asthmatics. In Cluster 2, the tendency for infected participants to be clustered together with their corresponding uninfected counterpart is statistically significant (p-value=0.007, 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, on average, the cophenetic distances (see Materials and methods section for the definition of cophenetic distance) between the infected cluster members and their uninfected counterparts are statistically significantly lower in Cluster 2 when compared to Cluster 3 (p-value=0.033, one-tailed Mann-Whitney-U-test, see Appendix 1-figure 3).
The sub-clusters found within Clusters 2 and 3, respectively (marked with orange and blue rectangles in Figure 1), were analyzed in terms of enrichment in or depletion of pre-and post-challenge time series. The results are presented in Table 3. This analysis provides evidence for a statistically significant separation of pre-and postchallenge 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 Appendix. The distance between any twotime 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 time series from asthmatics which are prominently different from those from other asthmatic subjects in Cluster three and also from healthy subjects in Cluster 2. These might be regarded as outliers.
Time series of percentage of eosinophils in nasal lavage fluid Findings are summarized in Table 4. The corresponding dendrogram is depicted in Appendix 1figure 1. In brief, three clusters were identified. Cluster 1 consists of four time series stemming from three asthmatics. As can be read off of the dendrogram depicted in Appendix 1-figure 1, 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, while the vast majority of the time series from healthy participants are grouped together in Cluster 2, most of the time series from asthmatic participants are split into two different Clusters, namely Clusters 1 and 3. This suggest a higher heterogeneity among the asthmatics. 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, on average, the cophenetic distances between the infected cluster members and their uninfected counterparts are statistically significantly lower in Cluster 2 when compared to Cluster 3 (p-value=8.96e-05, one-tailed Mann-Whitney-U-test, see Appendix 1-figure 2).

Autocorrelation properties of the biomarker time series
For every participant, the autocorrelation coefficient of the lung function parameters time series and of the FeNO time series was calculated using a for each parameter type physiologically meaningful time lag. More specifically, a one-day lag was used for lung function parameters, and a two-days lag for FeNO. Due to the low sampling frequency used before the viral challenge, the time series of eosinophil and neutrophil cell density in nasal lavage fluid were not included in the autocorrelation analysis.
The resulting autocorrelation coefficients were then used to compare the groups of asthmatic and healthy participants prior to and after the viral challenge (see Appendix 1-figures 10-14). Briefly, in terms of autocorrelation, the lung function parameters PEF (% predicted), FVC, and FEV1/FVC discriminate significantly between the pre-and post-infection status in both healthy and asthmatic participants. Indeed, for these three parameters, there is a moderate positive autocorrelation before the challenge, which then disappears after the challenge. Furthermore, in terms of autocorrelation, only FeNO discriminates significantly between healthy and asthmatic participants, and it does so only after the challenge. More specifically, only asthmatics exhibit a moderate autocorrelation of their FeNO time series after the challenge. However, after a multiple pairwise-comparison correction aimed at controlling the false discovery rate, only the parameters FVC and FEV1/FVC discriminate significantly between the pre-and post-infection status in both healthy and asthmatic participants, while the other findings lose their statistical significance. The results after multiple pairwise-comparison correction are summarized in Appendix 1-table 2.

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 of the biomarkers measured. However, prior to doing that, the efficacy of the inoculation with RV16 needed to be established. Indeed, 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 (Appendix 1-table 1). According to the above mentioned 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.
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, which regards time series as univariate empirical distributions of the biomarker at hand, aimed at detecting distributional changes in a given biomarker induced by the viral challenge: Here, each participant's pre-and post-challenge time series of each biomarker were compared using the Kolmogorov-Smirnov test. The second criterion aimed at detecting short-term and transient relative changes induced by the viral challenge in the context of the relative changes observed prior to the challenge. Here, throughout the entire period of observation, we assessed the relative change of each biomarker taking place within time intervals of 10 days. (see Subsection 5.2 and Figure 3 in the Materials and methods section below).

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 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. This loss of self-similarity is not merely the result of a larger response to the virus infection. Rather, we detected changes in the overall biomarker fluctuation dynamics elicited by the viral challenge that render asthmatics more similar to other infected asthmatic participants than to their uninfected counterparts.
Experimental evidence supporting our hypothesis of a loss of adaptive capacity of the respiratory system in asthma Our hypothesis in this study is based on the following 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 cluster analysis of the pre-and post-challenge time series of the percentage of eosinophils in nasal lavage fluid resulted in two main clusters: Cluster 2, which is statistically significantly enriched in healthy participants, and Cluster 3, which is mainly composed of asthmatic participants. In Cluster 2, the tendency for infected participants to be clustered together with their corresponding uninfected counterpart is statistically significant and clearly higher than in Cluster 3. In the clustering of the pre-and post-challenge 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. For both biomarkers (percentage of eosinophils in nasal lavage fluid and FeNO), the difference in this regard between Cluster 2 (mainly healthy participants) Table 3. Enrichment analysis of the sub-clusters found within the clusters described in Table 2 above (the clusters marked with orange and blue rectangles in Figure 1). Enrichment in pre-challenge time series is marked in bold letters, depletion of pre-challenge time series (and consequently enrichment in post-challenge time series) in italics; the corresponding p-values were calculated using the hypergeometric test. and Cluster 3 (mainly asthmatic participants) is further substantiated by the fact that, on average, the cophenetic distances between the infected cluster members and their uninfected counterparts are statistically significantly lower in Cluster 2 when compared to Cluster 3. In our clustering based on time series of cell density in nasal lavage fluid, there is a statistically significant separation of pre-and post-challenge time series. Moreover, within the subgroup of mainly post-challenge time series we found a cluster of size seven enriched in asthmatic participants (see Appendix 1- figure 4).
Summarizing, we have found evidence for the tendency of infected asthmatic participants to be more similar to other infected asthmatic participants than to their own uninfected counterparts when similarity is measured in terms of the biomarker dynamics of FeNO, of the percentage of eosinophils in nasal lavage fluid, and of cell density in nasal lavage fluid. This tendency was, however, not observed, when similarity was measured in terms of the biomarker dynamics of lung function parameters, and of neutrophil cell density in nasal lavage fluid.
Physiological interpretation of the group differences in autocorrelation properties of the biomarker time series We and others have previously shown that time series of lung function over days display a weak intrinsic autocorrelation over long-and short-time scales in asthmatic and healthy subjects (Frey et al., 2005;Delgado-Eckert et al., 2018;Thamrin et al., 2016). Lung function is correlated with lung function values of previous days; these correlation properties are related to severity of asthma and disease control (Thamrin et al., 2016). It has been hypothesized that such correlation properties are related to the balance between disease stability and adaptability of the system (Frey et al., 2011). We and others have previously shown that external stimuli such as medication can alter these correlation properties dependent on the applied drug action (Frey et al., 2005;Thamrin et al., 2016). Here, we provide first evidence that viral stimuli can also alter these correlation properties. In both healthy and asthmatic subjects, correlation of daily lung function was weakened by the viral challenge of the respiratory system. We hypothesize that following the viral challenge, the lung mechanical system properties are less deterministic consistent with a lower stability of the respiratory system. 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 Materials and 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. 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 Materials and methods section for the definition of neighbors. Asthma as a chronic disease: Causal chain of mechanisms or complex system behavior?
The classical analytical approach in asthma research is the identification of individual mechanisms (e. g., airway obstruction) or a series of mechanisms involved in the disease process (e.g. viral trigger fi inflammation fi bronchial hyperreactivityfi airway obstructionfi respiratory symptoms). However, epidemiological observations have questioned such a simple causal relationship between these mechanisms. For example, our previous work demonstrates that the strength of the trigger is often not proportional to the degree of response and symptoms (Frey and Suki, 2008). Also, the degree of inflammation, airway obstruction and bronchial responsiveness are often not closely related. Response to triggers (e.g., environmental pollutants) can occur with time lags and also in various degrees of intensity depending on the pre-existing conditions of the respiratory system (e.g. increased bronchial hyperresponsiveness following allergic sensitization or viral infection). Such a behavior is better explained by complex system behavior of a chronic disease. However, so far, the latter is difficult to capture and remained a theoretical concept.
We believe that the current experiment provided significant evidence to support the existence of such system behavior in asthma. In a well-established human challenge model, the respiratory system was challenged with a standardized (viral) stimulus, an established approach to test the behavior of complex network systems. We made the following observations: First, prior and after the challenge a set of lung functional and inflammatory asthma-biomarkers showed temporal fluctuations, in both healthy and asthmatic human subjects. There was a large inter and intra-individual variation. Second, the temporal relationship between inflammatory and lung functional biomarkers and symptoms was weak, not supporting the concept of simple proportional interactions of the above-mentioned causal chain of mechanisms. Third, despite these fluctuations and despite the absence of a strong response to the viral challenge in a distinct mechanistic biomarker, we identified differences between healthy and asthmatic humans in these dynamic variations. This indicates that the complex 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 postchallenge 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 10 day 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. interactions of inflammatory and lung functional parameters and thus the control of biological responses relevant to the respiratory system must be different in asthma. Our data support the hypothesis of a loss of adaptive capacity in asthma, which impedes the fast return to the pre-challenge stable dynamic steady state. We hypothesize that some clinical phenomena are consistent with such loss of adaptive capacity, such as, for example increased morbidity and prolonged respiratory symptoms after viral infection in asthmatics, persistent bronchial hyperresponsiveness after viral challenge, or slower return of airway obstruction following viral challenge (Busse et al., 2010).

Variable and heterogeneous effect of the viral challenge 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 preand 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 short-term/transient relative changes induced by the viral challenge that were statistically significantly different in magnitude from shortterm changes observed during the pre-challenge phase (see columns 4 and 5, and rows 1-4 in Table 5 above, and Supplementary 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 (Seemungal et al., 2001), 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 Table 6. The demographics of the study population. BMI is Body Mass Index. Only one healthy subject smoked two pack years or less 2 years before recruitment to our study, which is considered an insignificant smoking history. FEV1: forced expiratory volume in one second. PEF: peak expiratory flow.

Demographic features
Healthy Asthmatic nature (see rows 5-8 in Table 5 above, and Supplementary 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.
Potential physiological and inflammatory mechanisms responsible for the biomarker dynamics observed in the group of asthmatics and the resulting reduction in adaptive capacity One potential explanation for the differences in biomarker dynamics observed between healthy and asthmatic participants could be the level of airway obstruction. Indeed, Dames et al. quantitatively assessed the overall variability and complexity of airflow time series in patients with COPD during resting breathing. They found that airflow pattern complexity was reduced proportionally to airway obstruction measured with spirometric indices (Dames et al., 2014).
Changes at the level of airway smooth muscle (ASM) could also be associated with differences in the dynamics of system properties based on mechanical or environmental triggers affecting airway responsiveness (Brook, 2014;Noble et al., 2012). Lack of appropriate ASM-stretch has been associated with pathological incapacities of the asthmatic airways (Slats et al., 2007), which may mechanistically explain the inability to return back to homeokinetic range in asthma after a viral infection.
Inflammatory mediators and pathways involved in asthma could also explain the observed fluctuations in biomarker dynamics. With respect to the inflammatory biomarker dynamics it is well known that the production of inflammatory mediators is strictly controlled to ensure a limited but effective inflammatory response. To that end, most mRNAs encoding inflammatory mediators contain regulatory motifs like AU-rich elements in their 3'-untranslated region, which affect both transcriptional control, mRNA half-life and translational control (Hao and Baltimore, 2009). Changes in these regulatory mechanisms are likely to lead to less well-controlled inflammatory responses and in fact beyond that as the expression of many response genes are controlled by these regulatory motifs. We have recently shown that this translational control of AU-rich containing mRNAs in primary bronchial epithelial cells from mild and severe asthma patients compared to that from healthy controls is defective (Ravi et al., 2019). This led to exaggerated ex vivo production of, for example mediators driving neutrophilic inflammation, which correlated with the in vivo neutrophilic inflammation. In a follow-up study (A. Ravi et al, submitted) we showed that a RV16 challenge worsened this defect and with that the neutrophilic inflammation whereas this was not the case in healthy controls. The strong correlations between this defect and neutrophilic inflammation are suggestive of causality although formal proof is still lacking. This defect in the bronchial epithelial cells in asthma may underlie loss of adaptive capacity in response to a RV16 challenge.

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 and therefore our findings may not be directly translatable to moderate and severe asthma. A similar study setting for severe asthma is not feasible because it is ethically not acceptable to challenge those asthma patients with RV16 and in addition severe asthmatics are likely to be on corticosteroid treatment therefore introducing a confounding factor. Another possible 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.

Conclusion and implications
We presented 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 (Frey et al., 2011). As a proof of concept, such changes in homeokinetic system properties would provide evidence supporting the idea that not only singular factors in isolation, but their interaction and also system properties, such as the interactions between their constituent parts, may contribute to disease dynamics and phenotype stability. This systems-level understanding of chronic asthma may open up new avenues for better understanding of asthma and other chronic dynamic diseases. Already in this small sample size, it is obvious that there is remarkable individual temporal variability in inflammatory and physiological biomarkers, not only in disease but also in health. Thus, dynamic fluctuations of physiological processes around an equilibrium state, including their related biomarkers, are an intrinsic feature of the respiratory system. Moreover, even in health there are strong inter-individual differences in these dynamic characteristics. Nevertheless, within a given healthy subject, fluctuations remain similar following the virus challenge, indicating that, these dynamic fluctuations seem to be in a stable dynamic equilibrium state in health. This characteristic is lost in asthma.
Future studies involving time series of biomarker measurements may help us understand this system instability in chronic asthma or other airway diseases, 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/NL5317).

Participant cohort
Twelve non-smoking, atopic (as determined by positive skin prick test to common aeroallergens), 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 6. 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 7 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 Appendix). The study was mainly divided into two 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 two 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 Appendix 1-figure 15.

Measurement and collection of biomarkers 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 (Miller et al., 2005).
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.

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 (Grünberg et al., 2001) [Refer to Appendix for details]. For the calculation of short-term/transient changes, a gliding interval or window is moved, one day at a time, along the time series. The relative change between the first and last entry of the gliding window is calculated, resulting in a new time series of short-term relative changes r i . (B) A healthy participant's time series of short-term relative changes in FeNO is depicted. A gliding interval of size 10 days was used to calculate it from the participant's time series of FeNO measurements. 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 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. (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'. Table 8 provides an overview of the different sample measurements along with their frequency before and after rhinovirus challenge.

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 (Grünberg et al., 1999). All participants were screened for the presence of respiratory viruses just before the challenge, to rule out a concomitant infection resulting in a cold (see Appendix for more details). Those participants with a positive outcome of this test were excluded from the study. An experimental RV16 infection was induced by using a relatively low-dose inoculum of 100 TCID50 (Tissue Culture Infective Dose determining the amount of virus required to cause cytopathy in 50% of the cells) to mimic a natural exposure. The study protocol along with the viral dose used and its safety have been approved by the institutional Medical Ethics Committee in Amsterdam University Medical Centre, the details of which have been included in Appendix. Data from our previous study show that a low dose is sufficient to induce mild cold-symptoms (Mallia et al., 2006). 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 (Grünberg et al., 1999).
Refer to Appendix for further details.

Statistical and computational analysis
Statistical tests resulting in a p-value less or equal to 0.05 were regarded as significant.

Assessment of differences: Pre-vs. post-viral-challenge
For each participant, their time series of a given biomarker prior to and after the viral challenge were compared. This comparison was based on the Kolmogorov-Smirnov test, whereby the time series were treated as empirical distributions, thus disregarding the chronological order of the measurements. Differences in the variance between the pre-and post-challenge distributions were assessed using Levene's test (Levene, 1961).
Multiple comparison correction was performed where required, using the false discovery rate (FDR) method of Benjamini and Hochberg (1995), setting the expected proportion of falsely rejected null hypotheses to 0.05. Table 8. The overview of different measurements performed in the study along with the frequency of sampling before and after rhino-virus challenge. Measures 1-4 include repeated measurements and 5,6 represent one-time measurement to screen the subjects for the study. eight refers to the experimental intervention in the study. FEV1: forced expiratory volume in one second. FVC: forced vital capacity. PEF: peak expiratory flow. FeNO: fractional expired concentration of nitric oxide. The time series of a given biomarker, prior to and after the viral challenge, were regarded as empirical distributions and compared to each other using the Earth Mover's Distance (EMD) (Rubner et al., 1998). The resulting pair-wise distances between distributions were then used for hierarchical clustering of pre-and post-challenge distributions. See Figure 2 and the Appendix for more details. 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.

Measurements of biomarkers
Calculation of short-term/transient changes For each participant individually, and for each biomarker, throughout the entire period of observation, the biomarker's relative change in value taking place within time intervals of 10 days was calculated. This choice of time interval length was made based on published literature whereby 5 days post exposure to respiratory viruses was shown to be critical. Hence a 10 day window for comparison would include 5 days before challenge to contrast with 5 days after challenge (Denlinger et al., 2011). This was done throughout the entire period of observation considering all possible time intervals consisting of 10 consecutive days. In order to assess the statistical significance of the short-term relative changes possibly elicited by the viral challenge, the magnitude (that is, the absolute value) of relative changes observed during 10 day time intervals starting at least 10 days prior to the challenge were compared, by means of a Mann-Whitney-U-test, to the magnitude of relative changes that took place during 10 day time intervals that contained the day of the challenge. See Figure 3 and the Appendix for more details.

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. Statistically significant enrichment or depletion were established using the hypergeometric test (Rubner et al., 1998).
The relative location of leaves in the clustering dendrogram was quantitatively evaluated using the cophenetic distance (Sokal and Rohlf, 1962). 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 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 Appendix 1- figure 16).
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 above.
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 Figure 1 above). 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 Figure 1 above).

Autocorrelation of time series
The autocorrelation of a given time series was calculated using the sample Pearson correlation coefficient of the original time series and the time series resulting after forward-shifting the original time series by the lag utilized. In other words, if the original time series consists of the values x 1 ; . . . ; x n and L is the lag, all complete (i.e., no member of the pair is a missing value) pairs x t ; x tÀL ð Þ for t ¼ L þ 1; . . . ; n are used to calculate the covariance in the formula of the sample Pearson correlation coefficient. The funders had no role in study design, data collection and interpretation, or the decision to submit the work for publication. Appendix 1

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 one indicates the seroconversion of the antibodies against the virus at the end of the study, Column two reflects the response to the virus by the PCR product on the day three after challenge and finally Column three shows the symptoms developed in the volunteers after the viral inoculation). Positive response in either of the three categories was considered evidence for a successful viral challenge. Appendix 1-table 1. Effectiveness of the viral challenge in asthmatic and healthy participants, respectively. A participant is considered to have a successful viral inoculation if any one of the three tests is positive. one indicates positive response and 0 indicates a failed response in the corresponding tests indicated in columns.

Percentage of eosinophils in nasal lavage fluid
The corresponding dendrogram is depicted in Appendix 1-figure 1 below. tances between time series corresponding to infected healthy participants and their uninfected counterparts. Only time series belonging to Cluster two in the clustering dendrogram obtained using the percentage of eosinophils in nasal lavage fluid (see Appendix 1-figure 1 above) 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 three in the clustering dendrogram obtained using the percentage of eosinophils in nasal lavage fluid (see Appendix 1-figure 1 above) are contemplated here. The two distributions are statistically significantly different (p-value=8.96e-05, one-tailed Mann-Whitney-U-test, the cophenetic distances in Cluster three being, on average, higher than the ones in Cluster 2.

FeNO (Exhaled Nitric Oxide)
Appendix 1-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 two in the clustering dendrogram obtained using FeNo data (see Figure 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 three in the clustering dendrogram obtained using FeNo data (see Figure 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 three being, on average, higher than the ones in Cluster 2.

Cell density in nasal lavage fluid
We found one cluster of size 25 enriched in pre-challenge time series (p=0.041), and one cluster of size seven enriched in asthmatics (p=0.049). See Appendix 1- figure 2 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.
In 52.7% of the soft bootstrap dendrograms, the resulting empirical p-value was smaller or equal to 0.05.
Bootstrap distribution of the empirical p-values resulting from the permutation test used for establishing the statistical significance of the proportion of neighboring pre-and postpairs found in Cluster 3': Min. 1 st Qu. Median Mean 3rd Qu. Max. 0.000000 0.004287 0.014640 0.080993 0.080477 0.450480 In 57.5% of the soft bootstrap dendrograms, the resulting empirical p-value was smaller or equal to 0.05.
Soft bootstrapping of the time series of the percentage of eosinophils in nasal lavage fluid yielded comparable results (data not shown).

Autocorrelation analysis of time series of biomarkers
For every participant, the autocorrelation coefficient of the lung function parameters time series and of the FeNO time series was calculated using a for each parameter type physiologically meaningful time lag. More specifically, a one-day lag was used for lung function parameters, and a two-days lag for FeNO. Due to the low sampling frequency used before the challenge, the time series of eosinophil and neutrophil cell density in nasal lavage fluid were not included in the autocorrelation analysis.
The autocorrelation of a given time series was calculated using the sample Pearson correlation coefficient of the original time series and the time series resulting after forwardshifting the original time series by the lag utilized.
The figures below (S10-14) display the distribution of the autocorrelation coefficient for each biomarker compared (FEV1, FEV1/FVC, FVC, PEF and FeNO) in the two groups before and after the viral challenge. The table below summarizes the results of comparing, in terms of the autocorrelation coefficient, the groups of asthmatic and healthy participants prior to and after the viral challenge. This comparison was conducted using a two-way ANOVA, when the data fulfilled the statistical requirements for the conduction of such a test. The latter requirements are normality (tested for using the Shapiro test) and homoscedasticity (tested for using Levene's test). For those biomarkers or derived magnitudes that did not fulfill the condition of homoscedasticity, the groups were compared pair-wise using the t-test or the Mann-Whitney-Wilcoxon test, depending on whether normality and homoscedasticity were given or not. Moreover, the comparison between pre-and post-challenge groups was conducted using the paired versions of the aforementioned tests. When conducting all four pair-wise comparisons, p-values were corrected for multiple testing using the FDR method.
The tables below display each participant's individual pre-and post-challenge autocorrelation coefficient of the lung function parameters time series and of the FeNO time series. The accompanying p-values are the outcome of permutation tests aiming at assessing the statistical significance of the autocorrelation coefficients found. Patient IDs are indicated by Pxy, their health status using H/A, denoting Healthy or Asthmatic.
Appendix 1-table 3. Table displaying individual pre (0) and post (1) challenge autocorrelation coefficient at 1 day lag for FEV1 time series and associated p-values for the permutation tests conducted for assessment of the statistical significance of autocorrelation coefficients. Materials and Methods

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 five 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 FEV 1 !70% of predicted  . Airway hyperresponsiveness, indicated by a positive acetyl-ß-methylcholine bromide (MeBr) challenge with PC 20 9.8 mg/ml (Sterk et al., 1993) . 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 five 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 FEV 1 !80% of predicted Negative acetyl-ß-methylcholine bromide (MeBr) challenge or PC 20 !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 five 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. Conventionally, power calculations are based on comparison of groups using a single magnitude. This study was different in this regard, as multiple measurements were incorporated in a longitudinal setting to assess the temporal dynamics of different biomarkers before and after challenge with rhinovirus. Hence a proper a priori power estimation was not feasible due to not yet developed statistical tools. However, 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 (Turner et al., 2006a;Turner et al., 2006b). This sample size could probably miss capturing the changes in biomarkers with weak effect sizes. Nevertheless, this sample size 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 Appendix 1-figure 15.
Forced expiratory flow 25-75% (FEF 25-75) was measured. Patients were systematically instructed to perform the home monitoring PEF measurements. The electronic PEF devices were tested and compared to the lab spirometry measurements.

Home monitoring
The morning and evening lung function measurements at home was done by Peak Expiratory Flow Meter (PEF, Micro Diary, CareFusion).
The following routine lung function indices were measured along with Asthma Control Questionnaire FEV1: Forced Expiratory Volume in 1 st second FVC: Forced Vital Capacity FEV1/FVC: Ratio of the amount of air exhaled in 1 st second to the amount of air exhaled during a maximal expiration.
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-bmethylcholine bromide) according to the standardized tidal volume method that is operative in routine clinical diagnostics in the hospital (Sterk et al., 1993).

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 min with 50 ml PBS. 100 ml 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 ml and was loaded onto a single block.
The prepared slides were stained as soon as possible with Diff Quick one stain (Medion Diagnostics), dried and mounted on Depex.
3.1.6.2 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 [2000 g for 10 min (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.

Assessment of differences: Pre-vs. post-viral-challenge
Consecutive measurements of a given biomarker prior to and after the viral challenge resulted in pre-and post-challenge time series for each cohort participant. Except for the calculation of transient changes, the time series were treated as empirical distributions, thus disregarding the chronological order of the measurements.
The empirical distribution of any biomarker before and after viral challenge was compared using Kolmogorov Smirnov test. A participant is considered a responder with respect to a given biomarker if the outcome of the Kolmogorov-Smirnov test results in a p-value<=0.05. Differences in the variance between the pre-and post-challenge distributions were assessed using Levene's test.
The empirical distributions of a given biomarker, prior to and after the viral challenge, were compared to each other using the Earth Mover's Distance (EMD), see Figure 2 in the Main Manuscript. The resulting pair-wise distances between distributions were then used for hierarchical clustering of pre-and post-challenge distributions. Our clustering approach makes use of the entire collection (distribution) 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 subtler differences and similarities between the participants' measurements that are less likely to be captured by conventional methods based on averages. 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 (Gan et al., 2007) 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 (Rubner et al., 1998), and the agglomeration procedure was done according to Ward's minimum variance method (Gan et al., 2007). 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).

Calculation of short-term/transient changes
For each participant individually, throughout the entire period of observation, the relative change of each biomarker taking place within time intervals of 10 days was calculated. This was done throughout the entire period of observation considering all possible time intervals consisting of 10 consecutive days, see Figure 3 in the Main Manuscript. The rationale for this is as follows. Any changes taking place within a period of 10 days that contained the day of the viral challenge need to be interpreted and understood in the context of any changes taking place within a period of about 10 days during the entire pre-challenge phase of the study.
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., Bartels et al., 1995). These interpolated values were only used for the assessment of the short-term/transient response.
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 preor post-challenge distributions.
The relative location of leaves in the clustering dendrogram was quantitatively evaluated using the cophenetic distance (Sokal and Rohlf, 1962). 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 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. 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 Figure Appendix 1-figure 17 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.