PDF
1. Department of Mathematics and Statistics, Texas Tech University, Broadway and Boston, Lubbock, TX 79409-1042, USA.
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.
Background: A method comparison study is a topic of considerable interest in health and biomedical-related fields. The use of this study is to compare a new method with a standard established method. There are two typical steps in method comparison studies, namely modeling the data set and using the fitted model to analyze method comparison data.
Methods: As a usual practice to model method comparison data, many recommend a mixed-effects model which assumes constant error variance (homoscedasticity) and normality of error terms. However, these assumptions are generally violated in practice. Thus, in this study, our main goal is to propose a copula based model to deal with non-replicated method comparison data. Moreover, a simulation procedure is carried out to validate the proposed model by means of different copula models.
Results: Results indicate that the Normal and Gumbel copulas give more accurate results in terms of bias, mean-squared error and coverage probability values of estimators. Further it is confirmed that the accuracy of model increases with the Kendall tau (τ) correlation between methods and the number of observations. Furthermore, the proposed methodology is illustrated by analyzing finger and arm systolic blood pressure data. Besides the Total Deviation Index (TDI) and Concordance Correlation Coefficient (CCC) values are used to check the agreement between the two methods.
Conclusion: The proposed model based on the Copula method can be used to model the method comparison data with balanced and unbalanced data designs.
Keywords: Method comparison, Concordance Correlation, Copula theory, Mixed-effects model, Total Deviation Index
Quantifying the agreement between two methods of measuring continuous variables such as blood pressure, cholesterol level, heart rate, etc. is extensively used in health and biomedicalrelated disciplines. These methods may be a medical device, a clinical observer, an assay, etc. In general, new methods are compared with standard established methods. If two methods agree sufficiently well, we may prefer the cheaper or the less invasive method and, they can be used interchangeably. When introducing new methods, we consider whether they perfectly match with the standard or the most precise method, or whether it is possible to replace the new method with the old one. Method comparison studies are widely used to answer these types of questions. Therefore, the use of this study is to compare a new method with a standard established method [1,8,16].
Generally, there are two typical steps in the analysis of method comparison data. The first step is to model the data and the second step is to use the fitted model to perform inference on measures of agreement that indicate how close the methods are. Several methods have been proposed to model the method comparison data, among them the regression-based model [6,9] and the standard mixed-effects models [24] are the most important. A standard mixed-effects model is used predominantly to model method comparison data [2,5,25,26,30]. The validation of the mixed-effects model highly depends on assumptions such as constant error variance (homoscedasticity) and normality of error terms. However, these assumptions are violated in practice [2,5,30]. Basically, the error variability may change with the magnitude of measurements. But, dealing with this heteroscedasticity is important because the extent of agreement changes with the magnitude, and therefore is not constant. If heteroscedasticity is not considered, the resulting agreement evaluation will be misleading. Several models have been proposed in [25,26] to analyze the heteroscedastic measurements with replicates. But all these models do not accommodate non-replicated measurements. However, in practice there are several studies with non-replicated measurements, as taking replicated measurements is time consuming and expensive. Therefore, there is a need for a model to deal with heteroscedastic, non-replicated measurements. Thus, we propose a model using the copula theory where fewer assumptions are required to model non-replicated method comparison data.
In the modern world of statistics, the copula is widely used in quickly defining a correlation structure between two or more variables [20,21]. Copulas are multivariate distribution functions whose one-dimensional margins are uniform on the interval (0, 1). Simply, the copula is nothing more than a way of creating a joint probability distribution for two or more variables while preserving their marginal distributions [27]. From the definitions of copula, it is much clearer that the copula has proven to be a useful tool in the analysis of dependency structures in statistics [14,27]. Thus, the copula contains all information about the dependence between random variables [15]. Since most traditional representations of dependence are based on the linear correlation coefficient and restricted to multivariate elliptical distributions, the copula representations of dependence are free of such limitations [4]. In addition, the copula enables us to model marginal distributions and the dependence structure separately. Furthermore, most traditional measures of dependence are measures of pairwise dependence whereas copulas measure the dependence between all the random variables.
In this study, we propose a model with the use of the copula to deal with non-replicated method-comparison data. This article is organized as follows. In Section 2, we describe a copula-based model for heteroscedastic method comparison data. A simulation study to validate the proposed model is described in Section 3. Section 4 describes the agreement evaluation techniques based on the proposed model and Section 5 is used to show an application to blood pressure data using the proposed model. Finally, Section 6 concludes with a discussion. All the computations in this article have been performed using the statistical programming language R with using an inbuilt package called copula [22].
A copula is said to be a multivariate probability distribution for which the marginal probability distribution of each variable is uniform. The concept is first proclaimed by Sklar through Sklars theorem [27]. Let H be a joint distribution function with margins F and G, then there exists a copula C such that for all x, y in R,
If F and G are continuous, then C is unique; otherwise, C is uniquely determined on Range of F \ Range of G. Conversely, if C is a copula and F and G are distribution functions, then the function H defined by above is a joint distribution function with margins F and G. Further if U1,……,Un are U(0,1) then the function C from [0,1]… [0,1] into [0,1] defined by C(u1,……,un )=P{U1≤u1,………,Un≤un}.
More precisely our model can be defined as follows. Let X = (x1,x2,…,xn)T ndenotes the nx1 vector of measurements from method Y = (Y1,Y2,…,Yn)T, denotes the nx1 vector of measurements from method 2,
and are the predicted values of the methods 1 and 2 under the proposed copula model, E =(ε1, ε2,…,εn)T and Υ = (γ1,γ2,…,γn)T represents the error terms of each methods.
Where
and D1, D2 are the distributions of respectively and are the uniformly transformed variables of the distributions respectively. The is the copula parameter.
If the joint distribution is said to be some arbitrary distribution D then is defined as,
σc12 is the variance of the 1st method, σc22 is the variance of the 2nd method, μc1, μc2 represents the mean of the 1st method and mean of the 2nd method respectively ρc and is the Pearson correlation coefficient between the two methods.
One-parameter Archimedean copulas are listed in the Table 1. Archimedean copulas are popular because they allow modeling dependence in arbitrarily high dimensions with only one parameter, which governs the strength of dependence [3,10]. Since the copula mainly deals with extreme distributions having heavy tails, the traditional correlation coefficient (ρc) cannot describe the correlation structure between variables [17,23]. Therefore, rank correlation (τ) proposed by Kendall is used in this study. This is a non-parametric measure of association between random variables. This measure is mainly used to capture the association between ordinal data. The definition of this coefficient comes from concordant and discordant pairs of the random variables which are used to be checked the association [17]. The Kendall correlation coefficient has a simple closed form in terms of copula and is defined by [29] as
Table 1 : Parameter settings for the simulation study.
C(u,v) where is the copula of bivariate distribution of X and Y. For all elliptical copulas, especially Normal and Student t copulas, there is a relationship between linear correlation coefficient and Kendall tau correlation.
for Archimedean copulas Kendall tau correlation is related with copula parameter value. This Kendall’s tau correlation coefficient can be used for either interval or ordinal data [28]. This is more powerful than Pearson’s correlation, when data consist of outliers. Concordance measures assume their maximal (minimal) values if the concerning random variables are perfectly positive (negative) dependent. Dependence measures assume their extreme values if and only if the random variables are perfectly dependent. Moreover, these Concordance measures are zero in case of independence. Dependence measures are zero if and only if the random variables under consideration are independent.
This section comprises of evaluating and comparing parameters that are used in the proposed Copula based model. Several algorithms have been defined in [11-13] to carry out the simulation study by means of copula. These algorithms are based on Monte Carlo simulation and are applied to the five widely used single parameter copulas, four of which belong to the Archimedean family while the last one belongs to the elliptical family. The main interest throughout this simulation study is to check the accuracies of the point and interval estimators of copula parameter and the parameters used to perform the agreement evaluation of the proposed model in Section 2. To achieve this, we focused on the mean squared error(MSE), bias and coverage probability of these estimators.
This simulation procedure is carried out by considering different settings for copula parameter, number of observations (τ), and the 5 different copula types namely Frank(Fr), Clayton(Cl), Gumbel(Gu), Joe(Jo) and Normal(No). The combinations of these settings are used for the 1000 simulations. Table 2 briefly summarizes the parameter settings for the simulation study as below.
Table 2 : Parameter settings for the simulation study.
This simulation provides us with the estimated point and interval estimators of parameters of interest, which are tabulated in Table 3 that presents the estimated biases and MSE of the above estimators. The biases get higher value when the number of observations is low with value. Meanwhile, when number of observations gets larger, the biases becomes small in magnitude with positive and negative values. Most of time, the Normal copula contains the least bias values. It becomes further small with increasing (τ) value. In contrast, for small values, there is no any specific copula type which gives the minimum for all bias values. Although the values differ between parameters, Gumbel and Normal copula gives significantly minimum values. But the MSE values suggest the estimators under Normal copula and Gumbel copula gives superior values in terms of other values in the increased number of and , especially the MSE of copula parameter (θ) values.
Table 3 : Estimated biases and MSEs of estimators computed using five copulas. (Pa, Fr, Cl, Gu, Jo and No are Parameter, Frank, Clayton, Gumbel, Joe and Normal respectively).
Table 4 shows the estimated coverage probabilities of 95% confidence intervals computed by the five copulas. The entries for each type are almost identical with few significant differences. However, it is obvious that it is closer to 95% with large number of observations and correlation values. In contrast, the coverage probability σc1 of is around 80% for small n and θ.
Table 4 : Estimated coverage probabilities of 95% confidence intervals.
Concordance correlation coefficient and total deviation index are widely used agreement measures in method comparison studies [16].
Concordance Correlation Coefficient (CCC)
Concordance Correlation Coefficient (CCC) [18,19] is used to measure the agreement between two continuous measurements of methods. If the appropriate copula is c, then
where σc12 is the variance of the 1st method, σc22 is the variance of the 2nd method, μc1, μc2 represent the mean of the 1st method and mean of the 2nd method respectively and is the Pearson correlation coefficient between the two methods.
The concordance correlation coefficient ranges from -1 to +1. Perfect positive agreement is achieved when the coefficient value becomes +1. Meanwhile a perfect negative correlation is achieved when it becomes a negative 1. It is obvious that the concordance correlation coefficient becomes zero when these two methods have no agreement between them at all.
Total deviation index (TDI)
Total deviation index is defined as the p0th percentile of the absolute differences of assays if the appropriate copula is. Generally, 0.80≤p0≤ 0.95.
TDI values are non-negative and smaller values in TDI indicate better agreement, which becomes perfect when TDI is zero.
The blood pressure data set in [7] was used to illustrate the proposed method. The data set consist of 200 measurements, taken from two methods of measuring systolic blood pressure, namely a standard established method where measurements are taken by considering arm blood pressure, and a new method using finger blood pressure. Mercury column sphygmomanometer is used to measure both measurements. Blood pressure data set is modeled using the proposed copula based approach and the fitted model is then used to perform the method comparison studies. All these computations are performed by R programming language.
Stepwise model fitting
We used a stepwise procedure for model fitting. The initial step is the preliminary analysis where the behavior of the variables in the data set is studied. As the subsequent step, we fitted the copula based model. Next, we applied agreement evaluation techniques for making inference about two methods.
In Figure 1, the solid line shows the distribution of two variables arm blood pressure and finger blood pressure. AIC and BIC values were used to find the best distributions for the data set. Also, the AIC and BIC values of each fitted distributions is shown in Table 5.
Table 5 : AIC and BIC values of candidate distribution.
We wish to select the model that minimizes the information loss from the candidate models. By fitting several potential distributions and looking at their AIC and BIC values, it is concluded that log normal distributions fit the data best. Therefore, it was fitted to the original variables and the parameters were estimated by means of maximum likelihood estimation. Moreover, in Figure 1 we can figure out the original distribution and the fitted distribution. Figure 1 reveals that the log-normal distribution with dotted line curve fits better with the distribution of the original variables with solid line curve. The other curves do not give a closer estimation to the original distribution.
Figure 1 : The original distributions and the fitted distributions of the arm blood pressure and finger blood pressure measurements data respectively.
The P-P and Q-Q plots for each measurement were plotted for the further verification of the selected distribution. Figure 2 shows the P-P and Q-Q plots for arm blood pressure and finger blood pressure variables respectively. Since the entire scatter points are close to the reference line, these figures confirm the distribution that we have chosen is the best one. Figure 3 shows the values after probability integral transformation. The frequency plot in the left side shows the dispersion of the transformed variables, whereas the right side shows the corresponding P-P plot. From these plots, these transformed variables still follow the log normal distribution since it is reasonably well dispersed around the reference line. The transformed variables show a strong relationship between variables with Spearman correlation coefficient of 0.744. This accomplishes that the probability integral transformations have not changed the dependence structure of the variables.
Figure 2 : QQ plot of finger blood pressure measurements.
Figure 3 : (a) Histogram of arm blood pressure measurements.
Finally, the copula was used to form a bivariate distribution from these two variables. The parameter was estimated through the chosen maximum likelihood estimation procedure. Different types of copula such as Frank, Clayton, Gumbel, Joe, Normal and Student were used to create the bivariate distribution. To decide on the best copula, we used the AIC and BIC values, which were summarized in the Table 6.
Table 6 : AIC and BIC values for candidate copula.
The minimum AIC and BIC confirms the best copula, Gumbel copula with parameter 2.422 seems to fit the data well and hence, it was chosen to fit the bivariate distribution. To validate the fitted model, residual analysis was performed based on the estimated values and true values. It is easily seen that there is no any pattern but behave randomly. Thus, it suggests that the chosen model fits the data well.
Results of the application
Before the model fitting for the data for agreement evaluation, the variables were checked for the presence of heteroscedasticity. From the Bartlett’s test of variance (p-value=0.0361), it is confirmed that there is enough evidence to approve the heteroscedastic behavior of two methods Figure 4.
Figure 4 : Residual Plot of Copula based model.
Normal, Gamma, Log-normal, Weibull and Logistic distributions were fitted for each assay of measurements of finger and arm systolic blood pressure data set. Among the five distributions that have been fitted, the log normal distribution was selected as the best fit as it has the minimum AIC and BIC values. A similar step is carried out to select the best copula for the joint distribution function and the Gumbel copula is selected as the best fit for blood pressure data with parameter 2.42. To check the adequacy of the fitted model, the residuals were analyzed and results confirmed that this model fits the data well.
The detailed information on agreement measures is given in the Table 7. Since the CCC is 0.7913, there is a sufficient agreement between two methods of measuring systolic blood pressure using finger and arm methods. Moreover, the TDI of 25.46 ensures the better agreement between these two assays. Also, the 95% lower bound for CCC and upper bound for TDI were calculated. The lower bound for CCC is 0.7555, while the upper bound for TDI is 28.4911. Hence, the two methods agree sufficiently well to be used interchangeably. It is clearly noticeable that measuring blood pressure using finger sphygmomanometer is easier in practice. The difficulties of handling the arm systolic blood pressure are well mentioned in the paper [7]. Therefore, it is more convenient to use finger systolic blood pressure as an alternative to arm systolic blood pressure.
Table 7 : The estimates for CCC and TDI (with P0=0.90), their standard errors and 95% simultaneous confidence bounds (lower bounds for CCC and upper bounds for TDI).
We proposed a copula based methodology to deal with method comparison data. Results indicates that the normal copula gives more accurate results in terms of bias, MSE, and coverage probability values of estimators. It is also confirmed that accuracy of parameters increases with the Kendall tau cor-relation value between methods and number of observations. The proposed methodology is then used to analyze the blood pressure data. A good agreement between the two methods has been observed, using the estimated CCC and TDI values. In addition, this model does not assume the constant variance and normality of error terms. A limitation of our model is that it works only with non-replicated measurements. However, our methodology can be extended to work with replicated measurements.
CCC: Concordance correlation Coefficient
TDI: Total Deviation Index
The authors declare that they have no competing interests.
| Authors' contributions | KRA | LSN |
| Research concept and design | √ | √ |
| Collection and/or assembly of data | √ | √ |
| Data analysis and interpretation | √ | √ |
| Writing the article | √ | √ |
| Critical revision of the article | √ | √ |
| Final approval of article | √ | √ |
| Statistical analysis | √ | √ |
The authors are grateful to Dr. Alan Close for providing the systolic blood pressure data.
Editor: Feng Gao, Washington University School of Medicine, USA.
Received: 23-Apr-2017 Final Revised: 31-May-2017
Accepted: 05-Sep-2017 Published: 24-Sep-2017
Copyright © 2015 Herbert Publications Limited. All rights reserved.
