BioVar: an online biological variation analysis tool

Outlier detection

The first step of statistical evaluation of the BV data is outlier detection. This step should be carried out rigorously and thoroughly since even a single outlying observation can significantly influence the variability of the BV data components. It is recommended that this step should be performed using Cochran’s test and Reed’s criterion in three steps [3], [4]:

Normality

After the application of outlier detection, the normality assumption of the BV data should be checked. As recommended by Braga and Panteghini (2016), the normality control should be performed in two steps: (i) the normality of each measurand and (ii) the normality of mean values of all subjects [3]. The Shapiro–Wilk test is one of the most widely used and powerful tests for checking normality assumption [13]. The estimation of variance components using a nested ANOVA does not assume a normal distribution or any other distributions for the model effects [6], [14], [15]. However, calculation of confidence intervals and statistical significance tests usually require normality assumption [11].

Moreover, the direct calculation of RCV also requires normal distribution [16], [17]. An appropriate transformation, such as logarithmic or any other transformation method, can be applied to the BV data in the case of non-normal distribution. If a transformation applied to the BV data, the transformed-data should be re-tested for normality using the Shapiro–Wilk test.

Steady State

Since estimations of CV_I and RCV assume a steady-state condition [18], this assumption should be checked before performing statistical analysis for the BV data. This assumption requires that data from all subjects should be in a steady-state. In a steady-state condition, there should be no trend in the concentrations of measurand during the study period [16]. Linear regression analysis can be performed on mean [19] or median [20] group value over the whole study period to verify that all subjects are in the steady-state. Subjects are considered in the steady-state if the 95% confidence interval of the slope of the regression line contains 0 [19]. If this assumption does not hold, then a data transformation, such as multiples of the median (MoM) and its natural logarithm (lnMoM), can be used to establish a kind of steady-state condition [18]. The MoM is a simple measure that shows how far a single test result from the median of the sample results. Let assume Y_ijk is the value of k^th replicate for j^th sampling time from i^th subject, then, the MoM value for each measurement of Y_ijk can be calculated as follows:

Finally, natural logarithm of MoM (ln (MoM (Y_ijk))) is used to make the distribution of the pooled residuals Gaussian [16], [21].

Homogeneity

Another critical assumption is the homogeneity of within-subject and analytical variability (between replicates). The Bartlett [22] test can be used to check the homogeneity assumptions. The analysis should not be continued in the case of heterogeneity since the heterogeneity will cause estimations not to be generalizable to the whole population [11]. One alternative, in this case, is to remove outliers until homogeneity is obtained [20].

Subset analysis

The significant differences in the CV_I and CV_G between gender groups (i.e., female and male) must be controlled to determine whether the estimates should be reported separately for genders. If the 95% confidence intervals of CV_I and CV_G overlap between genders, then it is concluded that there are no significant differences between gender groups [19]. Therefore, all-subjects results should be reported for CV_I and CV_G. In the case of non-overlapping confidence intervals, the results of CV_I and CV_G should be reported separately for each gender.

Moreover, the differences between the mean values of genders can be compared using a Student’s t-test. If there is no significant difference between the mean values of males and females, all-subjects CV_I and all-subjects CV_G are used to calculate the analytical performance specifications. In the case of the significant difference between the mean values of genders, the lowest of the two CV_G estimate is used to calculate the analytical performance specifications [23].

Nested ANOVA

After detecting and removing outliers, controlling normality assumptions, checking the steady-state condition, testing homogeneity assumptions, and performing subset analysis, variance components (e.g., CV_G, CV_I, CV_A) can be estimated. A two-fold nested ANOVA is the most widely used type of method to estimate variance components for the BV data. Generally, a number of samples (j) for the BV data are collected from a certain number of subjects (i), and the samples are analyzed within replicates (k). A general structure of a nested ANOVA can be expressed as the following model:

where Y_ijk is the value of k^th replicate for j^th sample from i^th subject, μ is the population mean, G_i is the effect of CV_G, I_ij is the effect of CV_I, and A_ijk is the effect of CV_A. It is assumed that the variance components (G_i, I_ij, A_ijk) are mutually independent and uncorrelated random variables with a mean of zero and variances σ_G², σ_I², σ_A², respectively [16].

As previously mentioned in the Normality section, the calculation of confidence intervals and statistical significance tests as well as the direct calculation of RCV, require normality assumption. One option is to apply natural logarithmic (ln) transformation in the case of non-normal distribution. First, the raw data are transformed using ln transformation, and then the nested ANOVA method is applied to the transformed data. Finally, the transformed results are transformed back to the original scale.

Moreover, the estimation of CV_I is vital for diagnosing and monitoring patients and for setting analytical performance specifications [16]. A comprehensive simulation study performed by Roraas et al. (2016) showed that CV-ANOVA is the best-performing method, resulting in the least biased and narrowest range estimates based on the percentiles for the CVI compared to standard ANOVA and ln-ANOVA. In the CV-ANOVA method, each subject’s data is divided by that subject’s mean value, then a standard nested ANOVA procedure is applied [16] and this method can be applied independently of the data distribution [11], [16]. The CV-ANOVA, however, does not provide an estimate of CV_G, since each individual has a mean value of 1. Therefore, the CV_G can be estimated by the standard ANOVA [16].

Index of individuality and reference change values

The index of individuality (II) is a useful parameter that can be used to assess the utility of population-based reference intervals. The II compares the biological variation of a laboratory test within an individual to that between all individuals, and it can be calculated as follows [24]:

According to the study results by Harris (1981), if the II is <0.6, the measurand shows low individuality meaning that the values for an individual cover only a small part of the reference interval. Low II indicates the limited utility of population-based reference intervals. On the other hand, if II is >1.4, the measurand shows high individuality meaning that dispersion of values of each individual span much of the distribution of the reference interval [25]. Hence, population-based reference values are appropriate. When a measurand showing low II, it is recommended not to use conventional population-based reference intervals. Instead, it is suggested to use the reference change values (RCV) for monitoring of longitudinal time changes of measurand concentrations in an individual [3]. The RCV is described as the difference needed between two serial results from an individual that a change to be statistically significant [26]. The RCV can be calculated as follow:

where z_α is 1.96 for 0.05 significance level and z_α is 2.58 for 0.01 significance level. This method assumes a normal distribution and results in symmetrical limits. On the other hand, Fokkema et al. (2006) suggested calculation of asymmetric limits when model effects are assumed to be log-normal distributed [27]. Asymmetric limits for the RCV can be calculated as follow:

Here, CV_lnI is the CV on the log-transformed data, as CV_lnI = (ln (1+CV_I²))^1/2. Roraas et al. (2016) recommend the use of the asymmetrical limits regardless of the distribution of data [16].

Analytical performance specifications

There are two main sources of error for a measurand: random error and systematic error. Random error occurs in replicates and varies unpredictably, whereas systematic error remains constant or varies in a predictable manner [28]. The systematic error is evaluated by the bias, while the random error is assessed by the imprecision measured by the CV [29]. To measure the total error of a measurand, allowable bias and allowable imprecision are combined in a single metric as the total allowable error (TEa):

Here, 1.65 implies 95% of the results will within the TEa limit under normal distribution [30]. Fraser et al. (1997) proposed to calculate minimum, desirable, and optimal quality levels for allowable imprecision, allowable bias, and TEa as following [7]:

1. For imprecision:

   1. Optimum allowable imprecision <0.25 × CV_I

   2. Desirable allowable imprecision <0.50 × CV_I

   3. Minimum allowable imprecision <0.75 × CV_I

2. For bias:

   1. Optimum allowable bias <0.125 × (CV_I² + CV_G²)^1/2

   2. Desirable allowable bias <0.250 × (CV_I²+CV_G²)^1/2

   3. Minimum allowable bias <0.375 × (CV_I²+CV_G²)^1/2

3. For total error:

   1. Optimum TEa <1.65 × (0.25 × CV_I) + 0.125 × (CV_I² + CV_G²)^1/2

   2. Desirable TEa <1.65 × (0.50 × CV_I) + 0.250 × (CV_I² + CV_G²)^1/2

   3. Minimum TEa <1.65×(0.75 × CV_I) + 0.375 × (CV_I² + CV_G²)^1/2

The analytical performance specifications should be defined for each measurand to ensure that the measurement error does not affect the BV results. However, caution is needed when calculating the analytical performance specifications. As pointed out by Oosterhuis (2011), the maximum allowable bias and imprecision can be derived separately from the BV data, and to combine these two values into a single expression has no theoretical foundation and may lead to gross over-estimation of the TEa [30].

Web-tool development

A web-based tool is developed for the BV data analysis using the R language environment. This tool can be used to perform outlier detection, to control normality assumption, to check the steady-state condition, to test homogeneity assumptions, to perform subset analysis for genders, to carry out nested ANOVA, and to generate plots. The tool is developed using shiny package of the R language and other languages, including HTML, CSS, and jQuery used to polish the front-end. Several R packages were used both in the back-end and front-end of the tool. The tool has an interactive interface. Users can upload their datasets in .txt format. There is an available example dataset in the tool to help users to test the applicability of the tool. In an appropriate dataset, rows must represent subjects, while columns must represent variables, and the first row must include a header that indicates variable names. Analysis procedures are easy to perform as well. Users can select any desired analysis method, choose appropriate variables, set options, and execute the analysis from the sidebar panel of the tool. Desired outputs will be displayed on the main panel of the tool. Table and plot results will be displayed on separate sub-tabs on the main panel to avoid confusion. An automatic report can be downloaded in HTML file format. A help page can be found on the web-page of the tool. The web-tool can be accessed freely at https://turcosa.shinyapps.io/biovar/.

Dataset

We simulated a dataset based on R codes provided by Roraas (2017) [31] with a nested design following the scheme of 48 individuals (24 males and 24 females) with 20 samples (time points) measured in replicate. The CV_G, CV_I, and CV_A were set to 15, 5, and 1%, respectively. The simulated dataset uploaded to the tool as in a long format using the Data upload tab (Figure 1).

Figure 1:

Data upload tab of the tool. The dataset uploaded to the tool as in a long format.

After uploading the data set, we can now move on to the Analysis tab to perform the BV data analysis (Figure 2). There are separate sub-tabs in this section. Outlier detection results can be obtained using Cochran’s test and Reed’s criterion (Figure 2a). The normality assumption can be checked using the Shapiro–Wilk test (Figure 2b). The steady-state condition can be controlled by performing linear regression (Figure 2c). Homogeneity of analytical and within-subject variabilities can be tested using the Bartlett test (Figure 2d). Subset analysis can be performed by comparing 95% confidence intervals of CV_I and CV_G for males and females (Figure 2e). Nested ANOVA can be carried out to estimate variance components and analytical performance specifications for all-subjects and separately for males and females (Figure 2f). Mean and absolute range plots can be obtained after each outlier detection step (Figure 2g). An analysis report can be generated and downloaded in HTML file format (Figure 2h).

Figure 2:

Analysis tab of the tool. a) Outlier results can be obtained. b) The normality assumption can be controlled. c) The steady-state condition can be checked. d) Homogeneity assumptions can be tested. e) Subset analysis can be performed. f) Analysis of variance (ANOVA) components can be estimated. g) Mean and absolute range plots can be created. h) A report can be generated.

BV data analysis

First, outlier detection performed to detect outlying subjects using the Outliers tab. Three replicate measurements were detected as outliers using Cochran’s test (Table 1), and no subjects found as outlying subjects using Cochran’s test and Reed’s criterion. Accordingly, the three outlying replicates were removed from the dataset.

Next, we move on to checking the normality assumption using the Normality tab. First, we performed normality tests on sets of results from each individual using the Shapiro–Wilk test and found that all subjects follow normal distributions (all p values > 0.05). Then, we move on to check the normality assumption on mean values of the subjects using the Shapiro–Wilk test, and it is concluded that the data is assumed to be normal on mean values of subjects (p=0.998, Table 2)

Next, the Steady-State tab is used to check whether all subjects were in the steady-state condition during the study period. Linear regression analysis performed on the mean group values over the study period (Table 3), and the steady-state assumption was concluded since 95% confidence interval of the regression coefficient included 0 (95% CI: −0.015 – 0.014).

The Homogeneity tab is used to test the homogeneity of within-subject and analytical variability (between replicates). The Bartlett test is used to test the homogeneity assumptions, and the test results showed that the homogeneity assumptions are valid (Table 4) for both analytical variability (p=0.940) and within-subject variability (p=0.083).

The Subset tab is used to determine whether CV_I and CV_G results should be reported separately for gender groups. Since the 95% confidence intervals of both CV_I and CV_G overlap between males and females, it is concluded that all-subjects estimates should be reported for CV_I and CV_G (Table 5).

Next, the ANOVA tab is used to display (i) mean of the measurand with its confidence interval, II and the number of samples (Table 6), (ii) ANOVA table (Tables 7), (iii) analytical performance specifications (Tables 8) and (iv) RCV result with lower and upper limits (Table 9). For all subjects; CV_A is founded as 0.987% (95% CI: 0.928% – 1.054%), CV_I is 4.225% (95% CI: 3.952 – 4.536%), CV_G is 13.937% (95% CI: 11.577 – 17.487%), II is 0.303 and RCV is founded as 12.026% (11.292 – 12.866%).

The Plot tab can be used to create mean and absolute range plots (Figure 3). These plots are especially useful to inspect outlying subjects and compare subgroups (i.e., genders). Four plots will be displayed in this tab: (i) before outlier detection, (ii) after removing outliers among replicates using Cochran’s test, (iii) after removing outlier subjects using Cochran’s test, and (iv) after removing outlier subjects using Reed’s criterion. Finally, one can download the whole report using the Report tab. The report can be downloaded as an HTML format. The report includes all the results and the tables generated by the tool.

Figure 3:

Mean and absolute range plot after removing outliers among replicates using Cochran’s test.