Application of Probabilistic Multiple-Bias Analyses to a Cohort- and a Case-Control Study on the Association between Pandemrix™and Narcolepsy

Background An increase in narcolepsy cases was observed in Finland and Sweden towards the end of the 2009 H1N1 influenza pandemic. Preliminary observational studies suggested a temporal link with the pandemic influenza vaccine Pandemrix™, leading to a number of additional studies across Europe. Given the public health urgency, these studies used readily available retrospective data from various sources. The potential for bias in such settings was generally acknowledged. Although generally advocated by key opinion leaders and international health authorities, no systematic quantitative assessment of the potential joint impact of biases was undertaken in any of these studies. Methods We applied bias-level multiple-bias analyses to two of the published narcolepsy studies: a pediatric cohort study from Finland and a case-control study from France. In particular, we developed Monte Carlo simulation models to evaluate a potential cascade of biases, including confounding by age, by indication and by natural H1N1 infection, selection bias, disease- and exposure misclassification. All bias parameters were evidence-based to the extent possible. Results Given the assumptions used for confounding, selection bias and misclassification, the Finnish rate ratio of 13.78 (95% CI: 5.72–28.11) reduced to a median value of 6.06 (2.5th- 97.5th percentile: 2.49–15.1) and the French odds ratio of 5.43 (95% CI: 2.6–10.08) to 1.85 (2.5th—97.5th percentile: 0.85–4.08). Conclusion We illustrate multiple-bias analyses using two studies on the Pandemrix™-narcolepsy association and advocate their use to better understand the robustness of study findings. Based on our multiple-bias models, the observed Pandemrix™-narcolepsy association consistently persists in the Finnish study. For the French study, the results of our multiple-bias models were inconclusive.


Introduction
In June 2009 the World Health Organization (WHO) declared an influenza pandemic phase 6 [1] caused by a novel influenza A/H1N1 strain leading the WHO's Strategic Advisory Group of Experts recommending commercial-scale production of pandemic A(H1N1)pdm09 vaccines [2]. One of the vaccines authorized by the European Medicines Agency was Pandemrix™, an ASO3-adjuvanted, inactivated split-virion influenza virus A/California/7/2009 (H1N1)v-like strain vaccine (GSK Vaccines, Rixensart, Belgium). Worldwide, the vaccine was distributed in at least 47 countries, with over 30 million vaccinated individuals in Europe [3].
Early in 2010, an increase in narcolepsy cases was observed in children and adolescents in Finland and Sweden [4,5], with preliminary investigations suggesting a temporal association with Pandemrix™, the only pandemic influenza vaccine used in these countries. Narcolepsy is a rare chronic sleep disorder, typically characterized by episodes of sudden loss of muscle tone, sometimes triggered by strong emotions (cataplexy) [6]. Narcolepsy etiology appears to be linked to genetic (HLA-DQB1 Ã 0602 allele [7]) and environmental (e.g., streptococcal infection [8], H1N1 natural infection [9,10]) factors. To further assess the signal in Finland and Sweden, numerous observational studies at country level [11][12][13][14][15][16][17][18][19] were subsequently conducted in Europe and this in the best way possible under challenging circumstances (Table 1). Given the urgency of assessing the signal and the fact that the pandemic vaccine was no longer in use, all studies used retrospective designs combining data from various sources, thus potentially introducing systematic errors. The main sources of systematic errors that have been identified for the narcolepsy studies are ascertainment bias, disease onset misclassification, selection bias due to preferentially including exposed cases, diagnostic bias and confounding by indication and by natural H1N1 infection [20,21]. Additionally, studies were affected by limited power as a result of the rarity of narcolepsy disease.
The use of quantitative bias analyses to study the impact of such systematic errors is advocated by key-opinion leaders [22][23][24] and international health authorities including regulatory agencies [25] [26], especially for studies with far-reaching public health implications and leading to policy recommendations [22,23]. This was not systematically done in any of the narcolepsy studies. In most studies, the potential bias assessment was limited to a qualitative discussion.
Several methods exist to quantify the impact of uncontrolled confounding [22,24,[27][28][29][30][31] and bias [22,24,25,27] on the risk estimates. These methods involve bias parameters, for which the values can be chosen to be either wide-ranging and hypothetical ('Which bias levels are needed to reach a certain value?': target-adjustment sensitivity analysis) or evidence-based ('Which result would be obtained given a certain realistic level of bias?': bias-level sensitivity analysis) [29]. Traditional deterministic bias analyses estimate the adjusted effect measure in light of the observed data and a finite set of bias parameter values. Such analyses can be improved upon using (probabilistic) Monte Carlo simulation, randomly sampling from the probability distributions of the bias parameters to obtain a distribution of the adjusted effect measure [31]. To jointly evaluate a cascade of different bias sources, multiple-bias analyses has been proposed [22,24,25,27].
In this paper, we apply bias-level probabilistic multiple-bias analyses to two of the published narcolepsy studies in order to assess the robustness of the Pandemrix™ -narcolepsy association.

Selected studies
Multiple-bias analyses [22,24,25,27] provide bias-adjusted estimates for effect measures derived from standard 2x2 tables. From the studies providing such effect measures [12-14, 16, 17, 19], we initially selected the Finnish cohort study because it has the strongest risk estimate (i.e. the highest lower limit of the 95% CI). To illustrate multiple-bias analyses for other designs, we also selected the French study, being the only case-control study (Table 1).
Finnish pediatric cohort study. This population-based study compared narcolepsy rates in a cohort of vaccinated to unvaccinated children aged 4-19 years [16] between 1999 and 2010. Vaccination status was obtained from primary health care databases. Narcolepsy cases (n = 67) were ascertained from the national care register and registers of three specialized health care centers following the Brighton Collaboration (BC) criteria [32], levels 1 to 3. The cohort vaccination coverage was 75%. For the primary analysis, the index date was the first documented health care contact and the follow-up period ended on August 15, 2010, the day French case-control study. A multi-center case-control study was carried out in 14 French expert orphan disease narcolepsy centers [13]. 59 cases (all ages) were recruited from patients with narcolepsy-cataplexy (BC Levels 1 and 2) [32], between 1 October 2009 and 30 April 2011. 135 controls were recruited from: (i) patients from the hospitals to which the participating sleep centers belonged; and (ii) healthy volunteers from a national database (Narcobank). Vaccination status was obtained through telephone interviews. The total population vaccination coverage was 8.8%. H1N1 vaccines used in France were Pandemrix™ (72%) and Panenza™ (28%). For the primary analysis, the index date was the date of polysomnography-MLST diagnosis. When considering all ages, the odds ratio (OR) for Pandemrix™ was 5.5 (95%CI: 2.5-12.0).

Bias-level probabilistic multiple-bias analysis
For both studies, we built a Monte Carlo simulation model. The model structure consists of a sequence of bias correction equations [27] in the reverse order of the bias cascade [27,29]. The model subsequently accounts for exposure misclassification, disease misclassification, selection bias and uncontrolled confounding by age, risk group and natural H1N1 infection (model structures in S1 Table) with evidence-based parameters. Tables 2 and 3 give an overview of the input distributions for the Finnish and French study, respectively. S2 and S3 Tables give details on the derivation of the input parameters. For every analysis, we generated 500.000 runs. All analyses were carried out in R 3.0.1 [33] (example code in S4 Table).

Input parameters: Finnish pediatric cohort study
Exposure misclassification. We judged that recall bias was negligible as the index date was the first reported health care contact. The authors assessed incompleteness of the exposure data through reviewing the local vaccination records of 1000 randomly selected individuals from the study cohort [16]. This review revealed four discrepancies, all of whom had been vaccinated but were not recorded as such, implying an exposure sensitivity for non-cases (Se XjD = 0 ) of 0.995 given the cohort vaccination coverage of 75% [16]. To account for parameter uncertainty, Se XjD = 0 was assumed to follow a betapert distribution [34] A betapert distribution is a continuous, smooth, unimodal distribution described by a minimum (min.), most likely (mlik.) and maximum (max.) value and which is commonly used in risk analysis to reflect expert opinion or in situations with limited data. As the review did not reveal false positives, the exposure specificity for non-cases was assumed perfect (Sp XjD = 0 = 1). The review of the vaccination records of all narcolepsy cases did not show any discrepancies, implying perfect exposure sensitivity and specificity among cases (Se XjD = 1 = Sp XjD = 1 = 1).
Disease misclassification. Differential ascertainment and lack of blinding are two potential sources of disease misclassification. Differential ascertainment bias occurs if vaccinated cases are (1) more likely to seek health care, (2) more likely to be referred to a specialized center and/or (3) more likely to be diagnosed compared to unvaccinated cases. These factors act multiplicatively, implying large potential for differential ascertainment bias. Indirect evidence for differential ascertainment can be found in the Swedish studies, showing a decreasing risk estimate from 4.19 (95% CI: 1.76, 12.1) to 2.92 (95% CI: 1.78, 4.79) with one additional year of follow-up [18,35], likely explained by the increasing narcolepsy incidence in the unvaccinated cohort with the extended follow-up [18]. Additional evidence is the shorter time-to-diagnosis among vaccinated (12 months, 7.9 months, 6.9 months) compared to unvaccinated subjects (60 months, 47.6 months, 12.6 months) as observed in the French [13,17,19]. We obtained the disease sensitivity among the vaccinated (Se DjX = 1 ) and unvaccinated (Se DjX = 0 ) starting from the observed time-to-diagnosis intervals and assuming an exponential distribution for time. We further assumed that the disease sensitivities are highly correlated with a correlation of 0.95 and used Gaussian copula's to simulate correlated values. Bp(min, mlik, max) = betapert distribution with a minimum, most likely, and maximum value; I(c) = deterministic distribution, a probability distribution of a random variable which only takes a single value c; N(μ, σ2) = normal distribution with mean μ and variance σ2; RR CD = rate ratio between confounder and disease.
doi:10.1371/journal.pone.0149289.t002    OR CD = odds ratio between confounder and disease. + Because only disproportionality in sampling weights causes selection bias, the constant C t can take any value. doi:10.1371/journal.pone.0149289.t003 Lack of blinding among the referred potential cases may result in a higher false positive probability among vaccinated compared to unvaccinated subjects. In the Finnish cohort study [16], narcolepsy cases were classified with increasing uncertainty as BC level 1 (16%), level 2 (76%) and level 3 (8%). We assume that false positives occur among level 3 cases only. The differential false positivity rate due to lack of blinding was reflected by assuming that the proportion of false positives among the BC level 3 cases was higher in vaccinated (min. = 5%, mlik. = 35%, max. = 50%) compared to unvaccinated subjects (min. = 5%, mlik. = 15%, max. = 50%). We also assumed a correlation of 0.95 for the false positive rates.
Selection bias. Selection bias is negligible in this population-based cohort study because the Finnish computerized national register was used. In this register, all residents are recorded and assigned a personal identity code that remains unchanged throughout a person's lifetime.
Uncontrolled confounding. Confounding by age, by indication and by natural H1N1 infection are the most important recognized uncontrolled confounders. To correct for uncontrolled confounding, assumptions are needed about (i) the association between the confounder and the disease (RR CD ), (ii) the prevalence of the confounder among the vaccinated (P CjX = 1 ) and (iii) among the unvaccinated subjects (P CjX = 0 ) [24].
Age might have acted as a confounder because the typical narcolepsy onset is between 10-30 years of life and the vaccination coverage was different by age (5-9 years: 80.6%; 10-14 years: 81,9%15-19 years: 56,6%) [16]. For the bias analyses, we used the age groups 5-14 years and 15-19 years and derived the confounder-disease association (RR CD ) from the age-specific narcolepsy rates for the period 2000-2010 as obtained from a large European study on the incidence of narcolepsy [36]. The prevalence of the confounder (i.e. the age group 15-19 years) within the vaccinated (P CjX = 1 ) and unvaccinated Finnish population (P CjX = 0 ) was derived from the age-specific Pandemrix vaccination coverage [16].
Confounding by indication would occur if belonging to a particular risk group for which H1N1 vaccination was recommended also carried an increased risk to develop narcolepsy. We derived the confounder-disease association (RR CD ) from the British narcolepsy study [15], which reported the odds ratio with and without matching on risk group. Then, we estimated the distribution of 'risk group' within the vaccinated (P CjX = 1 ) and unvaccinated Finnish population (P CjX = 0 ) by assuming the vaccination coverage within the risk group being minimally the overall vaccination coverage of 75%.
Confounding by natural H1N1 infection would occur if individuals exposed to natural H1N1 infection prior to protection through vaccination would be more likely to get vaccinated (being more likely to seek care and be offered vaccination) compared to non-exposed individuals. We derived the confounder-disease association from a Chinese study, which showed strong seasonality in narcolepsy occurrences correlated with both H1N1 and seasonal influenza counts assuming a partial population exposure of 31.8% (95%CI: 29.1-34.1) as observed in a Chinese serological study among children (6-15yrs) [37]. We then calculated the prevalence of H1N1 infection among the vaccinated (P CjX = 1 ) and unvaccinated (P CjX = 0 ) Finnish population given the Finnish vaccination coverage of 75% while assuming an H1N1 attack rate of 29% (95% CI: 15.8-41.9) [38] and a relative risk of vaccination given H1N1 exposure of minimum 1.0, most likely 1.2 and maximum 1.5. Hereby, we assume, as in [39], that vaccination effectively took place after the first H1N1 season as the H1N1 vaccination campaign for the healthy population in Finland started a week after the peak epidemic in the country [40].
For each of the three recognized sources of uncontrolled confounding, we assumed a bivariate correlation of 0.95 between the confounder prevalence rates within the vaccinated and unvaccinated population.

Input parameters: French case-control study
Exposure misclassification. The self-reported vaccination status used within the study might induce recall bias [13]. For all cases, vaccination was documented, implying perfect exposure specificity (Sp XjD = 1 ). We assumed perfect exposure sensitivity (Se XjD = 1 ) since medical records of cases were scrutinized. For six controls, vaccination could not be documented [13]. We calculated the exposure specificity (Sp XjD = 0 ) assuming that all, three and none of the 'undocumented' controls were false positive. Exposure sensitivity among controls (Se XjD = 1 ) was assumed to be high, but not perfect. In studies comparing self-reported influenza vaccination status to vaccination records, the average exposure sensitivity is 96% [41]. Because a mass vaccination campaign was implemented in France during the H1N1 influenza period [42], the exposure sensitivity was assumed to be higher than 96%.
Disease misclassification. As only BC level 1 and 2 cases were retained, disease misclassification was assumed negligible.
Selection bias. Differential ascertainment (diagnosis is necessary for case selection) and participation bias potentially affect the selection probability of a case. Similar to the Finnish study and starting from the observed time-to-diagnosis intervals, we obtained the selection probability referring to the ascertainment of a vaccinated (P case(a)jX = 1 ) and unvaccinated case (P case(a)jX = 0 ). We obtained the selection probability due to participation for vaccinated (P case (b)jX = 1 ) and unvaccinated cases (P case(b)jX = 0 ) starting from the number of participating cases (n = 59) and the participation rate (71%). We assumed that exposure prevalence among the non-participating cases was minimum 30%, most likely 45% and maximum 59% (with 59% being the exposure prevalence among participating cases). As ascertainment bias and participation bias are independent given exposure, the overall sampling probability is the product of the bias-specific sampling probabilities. We further assumed that the selection probabilities related to the ascertainment of vaccinated and unvaccinated cases were highly correlated with a correlation of 0.95.
The selection probability of controls was judged to be subject to substantial uncertainty because (i) the participation rate among controls was not given, (ii) hospital controls were recruited from departments that specialized in treatment of patients for which H1N1-vaccination was not recommended and (iii) healthy controls were recruited from NARCOBANK, among which many healthcare workers for which vaccination was recommended [43]. Hence, the selection probability was unlikely to be independent of exposure with bias in either direction. Therefor, we assumed the ratio of the selection probabilities of vaccinated and unvaccinated controls (P controljX = 1 / P controljX = 0 ) to be minimum 0.8, most likely 1 and maximum 1.2.
Uncontrolled confounding. To adjust for confounding by age, we used the age groups 5-17 years, 18-29 years and 30-50 years having a H1N1 vaccination coverage of +/-10%, 4% and 7,5% [42]. Similarly as for the Finnish study, we derived the confounder-disease associations from the age-specific narcolepsy rates obtained from a large European study [36]. The prevalence rates of the polytomous confounder within the vaccinated and unvaccinated French population were derived from the age-specific H1N1 influenza vaccination coverage [16].
For confounding by indication, we made the same assumption for the confounder-disease association (OR CD ) as for the Finnish study and estimated the distribution of 'risk group' within the vaccinated (P CjX = 1 ) and unvaccinated French population (P CjX = 0 ) based on the French coverage data (risk group = 21%; non-risk group = 12%) [44].
Also for confounding by natural H1N1 infection, we made the similar assumptions as for the Finnish cohort study. We calculated the prevalence of H1N1 infection starting from the overall vaccination coverage in France (8.8%) and assuming a relative risk of vaccination given H1N1 exposure of minimum 1, most likely 1.2 and maximum 1.5.
Again, for each of the three recognized sources of uncontrolled confounding, we assumed a bivariate correlation of 0.95 between the confounder prevalence rates within the vaccinated and unvaccinated population.

Finnish pediatric cohort study
Given the assumptions made, we found disease misclassification (as a result of differential ascertainment) to have the largest decreasing impact on the rate ratio and confounding by age to have the largest increasing impact (Tables 4 and 5). Combining all recognized sources of bias and uncontrolled confounding, the original rate ratio of 13.78 (95% CI: 5.72-28.11) reduced to a median value of 6.06 (2.5 th -97.5 th percentile: 2.49-15.1). French case-control study In this study, we found selection bias (as a result of differential ascertainment, participation bias and improper selection of controls) to have the largest impact on the odds ratio (Tables 4  and 5

Discussion
We performed multiple-bias analyses on the Finnish [16] and French [13] narcolepsy studies using Monte Carlo simulation models with evidence-based distributions for the bias parameters. We are not aware that such systematic bias analyses embedded within formal models of bias and confounding were performed in any of the Pandemrix™-narcolepsy studies. Most sensitivity analyses reported so far assessed how certain operational study choices affected the results, such as using different onset dates [12,13,15,16], different follow-up periods [12,13,15,16] or different inclusion criteria [13,15]. Only in the UK study [15], deterministic singlebias analyses on case selection and confounder misclassification were carried out. In the Irish study [12], one deterministic bias analysis on case selection bias was performed. Although the risk estimates substantially decreased, they remained significant in these single-bias analyses.
In our multiple-bias analyses and given our base case assumptions, the Pandemrix™-narcolepsy association consistently persists in the Finnish study; The rate ratio of the Finnish study reduced from 12.7 to 6.06 (2.5 th -97.5 th percentile: 2.49-15.1) after adjusting for the recognized bias sources. The risk estimates of the French study proved to be less robust; Given our base case assumptions for all recognized bias source, the odds ratio decreased to 1.85 (2.5 th -97.5 th percentile: 0.85-4.08), with the inter-percentile range including 1. Ascertainment bias (i.e. disease misclassification in a population-based cohort and selection bias in a case-control study) was found to have the largest decreasing effect on the risk estimates. This large impact might be explained by the characteristics of narcolepsy, being rare under diagnosed disorder with often an insidious onset characterized by a very long interval from first symptoms to final diagnosis (up to 15 years on average [45]), in combination with a relatively short study period.
Our approach has limitations. We only quantified the impact of the recognized sources of bias and confounding, though other, so far unrecognized sources might exist. The assumptions were evidence-based to the extent possible. Nevertheless some assumptions were based on limited information, in particular the assumption on the H1N1 infection-narcolepsy association or the differential H1N1 attack rate among vaccinated and unvaccinated. We further assumed that the shortened time-to-diagnosis among the vaccinated cases [13,17,19] was an indication of ascertainment bias. This difference could alternatively be explained by a more severe clinical presentation. However, the comparison of other disease characteristics, such as hypocretin levels or sleep latency test results, between vaccinated and unvaccinated cases does not support this [17,19,43,46,47]. Another explanation for the differential time to diagnosis in [13,17,19] could be the differential observation periods, whereby unvaccinated cases could have symptom onset before the pandemic, but vaccinated cases could only have symptom onset started as of the pandemic. To assess this, we used an alternative parametrization of ascertainment bias, based on the Swedish observation that the narcolepsy incidence among unvaccinated more than doubled after public awareness [18]. This approach resulted in smaller reduction of the risk estimates with the association now consistently persisting for the Finnish and French study, with a median RR of 7.48 (2.5 th -97.5 th percentiles: 3.3-16.9) for the Finnish study and median OR of 2.26 (2.5 th -97.5 th percentiles: 1.09-4.67) for the French study when adjusting for all considered sources of bias.
Finally, it needs to be stressed that multiple-bias analyses assess the sensitivity of the observed risk estimate to hypothetical values of plausible sources of bias. They do not assess the presence of a specific source of bias, nor do they provide a true estimate of the causal association. Keeping in mind that conventional methods rely on the implicit assumption of no bias and no confounding, bias-analysis are a good quantitative extension of the qualitative speculations on bias and confounding that characterize a good discussion of study results [48]. We provided details on the derivation of the input parameters (S2 and S3 Tables) and example code (see S4 Table). This way, the analyses are fully transparent and reproducible and can easily be rerun using different distributions for the bias parameters.
We performed our analyses on two of the published Pandemrix™-narcolepsy studies and found that the association consistently persists in the Finnish study but not in the French study. We encourage the authors of the narcolepsy studies to take potential measured confounders into account as much as possible in their original statistical analyses and to perform multiple-bias analyses for these sources of bias and confounding that cannot be adjusted for in the original analyses. Whether the results reported in this paper change the conclusions of the Pandemrix™-narcolepsy association is a matter of debate. At the very least, more detailed sensitivity analyses allow a better understanding of the sensitivity of each study's finding to such presumed biases, which may eventually lead to more informed benefit-risk and public health decision-making.

Trademark statement
Pandemrix is a trademark of GSK group of companies.
Supporting Information S1