Method comparison with repeated measurements -- Passing-Bablok regression for grouped data with errors in both variables

TL;DR

This paper extends the Passing-Bablok regression to grouped data with errors in both variables, providing an unbiased estimator and a method comparison test with higher power for repeated measurements.

Contribution

It introduces the Block-Passing-Bablok estimator for grouped data with errors in both variables and derives its asymptotic properties and variance.

Findings

The estimator is asymptotically normally distributed.

Variance depends on measurement overlap between groups.

The method improves power in method comparison with repeated measurements.

Abstract

The Passing-Bablok and Theil-Sen regression are closely related non-parametric methods to estimate the regression coefficients and build tests on the relationship between the dependent and independent variables. Both methods rely on the slopes of the connecting lines between pairwise measurements. While Theil and Sen assume no measurement errors in the independent variable, the method from Passing and Bablok accounts for errors in both variables. Here we consider the case where multiple, e.g. repeated, measurements with errors in both variables are available for m samples. We show that in this case the slopes between repeated measurements need to be excluded to obtain an unbiased estimate. We prove that the resulting Block-Passing-Bablok estimate for grouped data is asymptotically normally distributed. If measurements of the independent variable are without error the variance of the estimate equals the variance of the Theil-Sen method with tied ranks. If both variables are measured with imprecision the result depends on the fraction of measurements between groups that fall within the range of each other. Only if no overlap between measurements of different groups occurs the variance equals again the tied ranks version. Otherwise, the variance is smaller. We explicitly compute this variance and provide a method comparison test for data with repeated measurements based on the method from Passing and Bablok for independent measurements. If repeated measurements are considered this test has a higher power to detect the true relationship between two methods.

Tables

Table 1. Table 1: We consider four different combinations of group sizes denoted by 100-100 ( p 1 = p 2 = 100 ) subscript 𝑝 1 subscript 𝑝 2 100 (p_{1}=p_{2}=100) , 180-20 ( p 1 = 180 , p 2 = 20 ) formulae-sequence subscript 𝑝 1 180 subscript 𝑝 2 20 (p_{1}=180,p_{2}=20) , 10 × \times 100 ( p 1 = ⋯ = p 10 = 100 ) subscript 𝑝 1 ⋯ subscript 𝑝 10 100 (p_{1}=\dots=p_{10}=100) , and 820-9 × \times 20 ( p 1 = 820 , p 2 = ⋯ = p 10 = 20 ) formulae-sequence subscript 𝑝 1 820 subscript 𝑝 2 ⋯ subscript 𝑝 10 20 (p_{1}=820,p_{2}=\dots=p_{10}=20) . The true values of the samples are given by x ~ k = k subscript ~ 𝑥 𝑘 𝑘 \tilde{x}_{k}=k and y ~ k = b ​ x ~ k subscript ~ 𝑦 𝑘 𝑏 subscript ~ 𝑥 𝑘 \tilde{y}_{k}=b\tilde{x}_{k} , for k = 1 , … , m 𝑘 1 … 𝑚 k=1,\dots,m . In addition, we vary the variance within groups, setting σ ϵ = σ η subscript 𝜎 italic-ϵ subscript 𝜎 𝜂 \sigma_{\epsilon}=\sigma_{\eta} to either 0.2 0.2 0.2 or 0.4 0.4 0.4 to create high or low overlap between groups. We simulate the performance under the null hypothesis β = 1 𝛽 1 \beta=1 and alternative scenarios where β ∈ { 0.98 , 0.8 , 0.2 } 𝛽 0.98 0.8 0.2 \beta\in\{0.98,0.8,0.2\} . The difference between the classical Passing-Bablok regression (cPBR) and the Block-Passing-Bablok regression (BPBR) from Definiton  2.6 is shown. In particular, we provide estimations for 𝔼 ​ [ β ^ ] 𝔼 delimited-[] ^ 𝛽 \mathbb{E}[\hat{\beta}] , 𝔼 ​ [ I ] 𝔼 delimited-[] 𝐼 \mathbb{E}[I] , the probability that the true slope β ∈ I 𝛽 𝐼 \beta\in I and the probability to reject the null hypothesis β = 1 𝛽 1 \beta=1 .

slope	groups	overlap	mean $\hat{b}$		mean $I$		$\mathbb{P}\left(\beta \in I\right)$		$\mathbb{P}\left(1\notin I\right)$
			cPBR	BPBR	cPBR	BPBR	cPBR	BPBR	cPBR	BPBR
$\beta =1.0$	100-100	low	1.001	1.001	[.923,1.085]	[.923,1.085]	.950	.950	.050	.050
		high	1.003	1.004	[.870,1.156]	[.856,1.177]	.950	.949	.050	.051
	180-20	low	1.003	1.002	[.860,1.169]	[.874,1.148]	.951	.949	.049	.051
		high	1.005	1.01	[.812,1.246]	[.771,1.326]	.951	.947	.049	.053
	10$\times$100	low	1	1	[.994,1.006]	[.994,1.006]	.950	.950	.050	.050
		high	1	1	[.988,1.012]	[.987,1.012]	.952	.952	.048	.048
	820-9$\times$20	low	1	1	[.988,1.013]	[.991,1.009]	.951	.952	.049	.048
		high	1	1	[.977,1.024]	[.982,1.018]	.951	.951	.049	.049
$\beta =0.98$	100-100	low	.983	.981	[.907,1.067]	[.904,1.064]	.950	.950	.071	.076
		high	.989	.984	[.856,1.141]	[.838,1.156]	.950	.948	.054	.057
	180-20	low	.991	.983	[.850,1.157]	[.856,1.128]	.949	.947	.053	.061
		high	.998	.991	[.804,1.239]	[.753,1.306]	.948	.951	.049	.051
	10$\times$100	low	.980	.980	[.974,.986]	[.974,.986]	.950	.949	1	1
		high	.980	.980	[.968,.993]	[.968,.993]	.950	.951	.876	.880
	820-9$\times$20	low	.981	.980	[.969,.994]	[.971,.989]	.950	.951	.838	.994
		high	.982	.980	[.959,1.006]	[.963,.998]	.946	.948	.326	.607
$\beta =0.8$	100-100	low	.828	.801	[.757,.906]	[.730,.877]	.888	.953	.987	.998
		high	.854	.807	[.730,.997]	[.672,.963]	.881	.952	.536	.679
	180-20	low	.888	.802	[.754,1.051]	[.686,.934]	.769	.949	.303	.823
		high	.924	.812	[.738,1.158]	[.592,1.095]	.770	.950	.117	.305
	10$\times$100	low	.801	.800	[.795,.807]	[.794,.806]	.934	.948	1	1
		high	.802	.800	[.791,.813]	[.789,.812]	.935	.947	1	1
	820-9$\times$20	low	.813	.800	[.802,.825]	[.792,.808]	.398	.953	1	1
		high	.824	.800	[.803,.847]	[.784,.816]	.385	.948	1	1
$\beta =0.2$	100-100	low	.317	.202	[.257,.383]	[.144,.261]	.022	.948	1	1
		high	.461	.279	[.352,.588]	[.165,.404]	.001	.758	1	1
	180-20	low	.559	.201	[.415,.753]	[.105,.302]	0	.952	.975	1
		high	.508	.280	[.318,.733]	[.093,.501]	0	.906	.685	1
	10$\times$100	low	.204	.200	[.200,.209]	[.196,.205]	.567	.954	1	1
		high	.212	.204	[.203,.221]	[.195,.213]	.236	.846	1	1
	820-9$\times$20	low	.274	.200	[.255,.302]	[.194,.206]	0	.951	1	1
		high	.328	.202	[.299,.366]	[.189,.215]	0	.941	1	1

Equations

$${\tilde{y}}_i=\alpha +\beta {\tilde{x}}_i,$$

$$x_i={\tilde{x}}_i+{\epsilon }_i,$$

$$y_i={\tilde{y}}_i+{\eta }_i,$$

$${\tilde{y}}_k=\alpha +\beta {\tilde{x}}_k.$$

$$x_{i,k}={\tilde{x}}_k+{\epsilon }_{i,k},$$

$$y_{i,k}={\tilde{y}}_k+{\eta }_{i,k},$$

$$\underset{1\le i\le p_k}{\max }{x}_{i,k}<\underset{1\le j\le p_l}{\min }{x}_{j,l}\forall k,l\in \{1,\dots m\},k<l,$$

$$\mathbb{E}[S_{i,j}]=\beta (1-\mathbb{E}[\frac{{ϵ}_i-{ϵ}_j}{{\tilde{x}}_i-{\tilde{x}}_j+{ϵ}_i-{ϵ}_j}])$$

$$S_{ij}^{kl}$$

$$K=\#\{S_{ij}^{kl}|S_{ij}^{kl}<-1\}.$$

$$\displaystyle\hat{\beta}=\begin{cases}S_{(\frac{N+1}{2}+K)},&\text{if $n$ is odd},\\ \frac{1}{2}\big{(}S_{(\frac{N}{2}+K)}+S_{(\frac{N}{2}+K+1)}\big{)},&\text{if $n$ is even}\end{cases}$$

$$\hat{\alpha }=\text{median}\{y_i-\hat{\beta }{x}_i\}$$

$$S_{11}^{kl}=\frac{y_{1,k}-y_{1,l}}{x_{1,k}-x_{1,l}}\text{for}1\le k<l\le n.$$

$$P(\beta )=$$

$$Q(\beta )=$$

$$\tilde{C}:=P(\beta )-Q(\beta ).$$

$$S_{i,j}^{k,l}=\beta \frac{x_{i,k}-x_{j,l}+z_{i,j}^{k,l}}{x_{i,k}-x_{j,l}+z_{i,j}^{{}^{\prime }k,l}}$$

$$\mathbb{V}[\tilde{C}]=\frac{1}{18}\bigg{(}n(n-1)(2n+5)-\sum\limits_{k=1}^{m}p_{k}(p_{k}-1)((2p_{k}+5)+4\sum\limits_{u\neq k}^{m}p_{u}q_{ku})\bigg{)}.$$

$$\mathbb{V}[\tilde{C}]=\frac{n}{18}\big{(}3(n-p)+2(n^{2}-p^{2})\big{)}-\frac{2}{9}p^{2}(p-1)\sum\limits_{k=1}^{m}\sum\limits_{u\neq k}^{m}q_{ku}.$$

$$\mathbb{V}[C]=\frac{n(n-1)(2n+5)}{18}.$$

$$\frac{1}{18}\bigg{(}n(n-1)(2n+5)-mp(p-1)((2p+5)-4p^{2}(1-p)\sum\limits_{k=1}^{m}\sum\limits_{u\neq k}^{m}q_{ku})\bigg{)},$$

$$\mathbb{V}[\tilde{C}]\sim\frac{1}{9}n^{3}\Big{(}1-l_{m}-l_{o}\Big{)}\text{, with}$$

$$l_m$$

$$l_o$$

$$\mathbb{V}[\tilde{C}]\sim \frac{1}{9}{n}^3(1-\frac{1}{m^2}).$$

$$\tilde{C}_{\gamma}=w_{\frac{\gamma}{2}}\tilde{\sigma},\quad M_{1}=~{}\Big{\lfloor}\frac{N-\tilde{C}_{\gamma}}{2}\Big{\rfloor},\text{ and}\quad M_{2}=~{}N-M_{1}+1.$$

$$\displaystyle I=\big{[}S_{(M_{1}+K)},S_{(M_{2}+K)}\big{]}.$$

$$a_L=$$

$$a_U=$$

$$\mathbb{V}[\tilde{C}]=\frac{1}{18}\bigg{(}n(n-1)(2n+5)-\sum\limits_{k=1}^{m}p_{k}(p_{k}-1)(2p_{k}+5)\bigg{)}.$$

$$\mathbb{V}[\tilde{C}]=\frac{n}{18}\big{(}3(n-p)+2(n^{2}-p^{2})\big{)}.$$

$$\tilde{C}:=\#\{S_{ij}^{kl}|S_{ij}^{kl}>\beta \}-\#\{S_{ij}^{kl}|S_{ij}^{kl}<\beta \}$$

$$Z_{ij}^{kl}:=\{\begin{array}{ll}sgn(S_{ij}^{kl}),& k\ne l\\ 0,& k=l\end{array}$$

Full text

keywords: Passing-Bablok regression, Theil-Sen regression, method comparison, method agreement, non-parametric statistics, rank-statistics, non-continuous sampling distribution, repeated observations, errors-in-variable, clinical chemistry

1 Introduction

Consider the simple linear regression, where $y=\alpha+\beta x$, for $x,y\in\mathbb{R}$. Compared to ordinary least squares regression (OLS), non-parametric estimation of the regression coefficients $\alpha$ and $\beta$ can be more robust. Deviations from normal distributed error terms, like heavy tails or outliers, do not disturb the estimation as strongly as for OLS. A good introduction to robust estimation is given in [13].

Here we consider the non-parametric Theil-Sen regression (TSR) and Passing-Bablok regression (PBR), which are both based on Kendall’s rank correlation [4] and provide a robust estimate of $\beta$. PBR and TSR do not rely on the assumption of normally distributed errors. The Theil-Sen estimate [11, 10] for $\beta$ is given by the median of all slopes of the connecting lines between pairwise measurements. If measurement errors occur for both $x$ and $y$, the TSR is biased towards zero. This phenomenon is also known from OLS and called regression dilution or attenuation. For least squares, Deming regression instead of OLS can be used to account for errors in both variables. In addition, [2] suggest how to use repeated measurements to correct for regression dilution, if errors are normally distributed. A variation of the Theil-Sen regression that accounts for errors in both variables is the Passing-Bablok estimator [8, 9]. The Passing-Bablok estimate is also given by the median of all slopes, but shifted by an offset $K$ to ensure that $x$ and $y$ are interchangeable. Interestingly, PBR and TSR seem to be popular in separate fields. While TSR is popular in Metrology and Environmental Science, PBR is mainly used in Clinical Biochemistry, Pharmacology and Laboratory Medicine to compare two alternative measurement methods. This might be due to the fact that PBR is outlined in a guideline of the Clinical and Laboratory Standards Institute [7]. In [3] a protocol for method comparison studies in clinical laboratories is suggested, including a paragraph on PBR. Furthermore, there have been attempts to use PBR for batch effect removal in gene expression analysis [6]. Throughout this manuscript we will refer to the Passing-Bablok regression, but note that setting $K=0$ in the PBR equals the TSR and the presented results naturally hold for both methods.

We consider the case, where $p_{k}$ repeated measurements of two methods, for $m$ true pairs $(\tilde{x}_{k},\tilde{y}_{k})$ are available. The alternative measurements of the two methods are given by data points $(x_{i,k},y_{i,k})$, $i\in\{1,\dots,p_{k}\}$. This scenario includes simultaneous repeated measurements of the same sample with two methods, as well as multiple measurements of the same sample separately measured for each method. In the later, measurements need to be combined randomly into observation pairs $(x_{i,k},y_{i,k})$. If the data set for a PBR contains repeated measurements it is important to account for this. In contrast to the expected slope $\beta$ between independent measurements, the expected slope between repeated measurements is [math]. Hence, the slopes between points that correspond to repeated measurements of the same underlying true value are meaningless and would distort PB estimates and lower the power of the associated statistical test if included. Here we describe how the variance changes if we omit the meaningless slopes between repeated measurements and provide the resulting test on the equivalence of two methods.

The paper is organized as follows: First, we consider the assumptions and definitions of Passing and Bablok and introduce the necessary adaptation for repeated measurements in section 2. Section 3 states our results, including the asymptotic confidence interval for the estimated parameter $\hat{\beta}$ in Corollary 3.5. The test for the equivalence of two methods with repeated measurements is given in Remark 3.7. In section 4 we discuss the implications of the suggested Block-Passing-Bablok procedure and compare it empirically to the PBR without repeated measurements. Finally, proofs are given in section 5.

2 The Passing-Bablok regression for repeated measurements

2.1 The standard Passing-Bablok regression for independent measurements

Under the hypotheses of a structural linear relationship between two measurement methods, i.e.

[TABLE]

for $\tilde{x}_{i},\tilde{y}_{i}\in\mathbb{R}$, where $y$ corresponds to one and $x$ to the other method, Passing and Bablok [8] considered the following assumptions:

Assumptions 2.1**.**

(standard Passing-Bablok Regression)

• i)

All points $(x_{i},y_{i})$ are of the form

[TABLE]

with $\tilde{x}_{i},\tilde{y}_{i}$ the (non-random) ’true’ values of the measurements and error terms $\varepsilon_{i},\eta_{i}$ and $i\in\{1,\dots,n\}$.

• ii)

$\varepsilon_{1},\dots,\varepsilon_{n}$ are iid and come from an arbitrary continuous distribution with mean zero, and $\eta_{1},\dots,\eta_{n}$ are iid and come from the same type of distribution such that $\beta\varepsilon_{i}$ and $\eta_{i}$ are iid. Moreover, $\varepsilon_{i}$ and $\eta_{j}$ are independent for each $i,j$.

2.2 Block-Passing-Bablok regression for multiple dependent measurements

Sometimes the data at hand is not independent as stated in assumption ii)(ii) but grouped into dependent sets of measurements. Examples include common situations, e.g. if measurements have been repeated multiple times for the same sample or a series of measurements is done under the same conditions. In this case, the Passing-Bablok Regression cannot be used straight away. Here, we expand the Passing-Bablok Regression such that it considers the data to be available in $m$ groups with $p_{k}(k=1,\dots,m)$ members. We index our data points $(x_{i,k},y_{i,k})$ by group $k\in\{1,\dots m\}$ and individual $i\in\{1,\dots p_{k}\}$ from this group. Each $(x_{i,k},y_{i,k})$ represents a measurement of the true values $(\tilde{x}_{k},\tilde{y}_{k})$ and we assume again a linear relationship

[TABLE]

We make the following assumptions for repeated measurement data:

Assumptions 2.2**.**

(Passing-Bablok Regression for grouped data)

• i)

All $n=\sum_{k=1}^{m}p_{k}$ points $(x_{i,k},y_{i,k})$ are of the form

[TABLE]

with $\tilde{x}_{k},\tilde{y}_{k}$ the ’true’ values of the measurements in group $k$ and error terms $\varepsilon_{i,k},\eta_{i,k}$ indexed by group $k\in\{1,\dots m\}$ and individual $i\in\{1,\dots,p_{k}\}$ from this group.

• ii)

All $\varepsilon_{i,k}$ are iid and come from an arbitrary and continuous distribution with mean zero, and all $\eta_{i,k}$ are iid and come from the same type of distribution such that $\beta\varepsilon_{i,k}$ and $\eta_{i,k}$ are iid. Moreover, $\varepsilon_{i,k}$ and $\eta_{j,l}$ are independent for each $i,j,k,l$.

Assumptions 2.3**.**

(non overlapping groups)

All $m$ groups are strictly separated on the x-axis, i.e. almost surely

[TABLE]

Remark 2.4**.**

Assumption 2.3 can only be fulfilled if the error distribution is not supported on the whole real line, e.g. for uniformly distributed errors if $|\tilde{x}_{i}-\tilde{x}_{j}|$ is large enough for all $i\neq j$. Note further that groups need to be also separated at the $y$-axis as in Assumption 2.3 to maintain the symmetry of the method. However, we will mainly consider the general case with possibly overlapping groups such that only Assumptions 2.2 i)-ii) hold.

Remark 2.5**.**

Assumption ii) ii) and 2.2 ii) set the measurement error ratio equal to the regression slope parameter $\beta$, i.e. for the variances of error terms, $\frac{\sigma_{\eta}^{2}}{\sigma_{\varepsilon}^{2}}=\beta^{2}$ holds. The assumption ensures that the slopes are symmetrically distributed around $\beta$ (see remark 3.2). A common case where the assumptions are fulfilled is the comparison of two measurement methods under the null hypothesis that both methods are equivalent, i.e. $\beta=1$, and all $\epsilon_{i}$ and $\eta_{i}$ are iid.

Estimating the regression parameter $\beta$ from grouped data

The Passing-Bablok regression for independent measurements [8] makes use of the slopes $S_{i,j}=\frac{y_{i}-y_{j}}{x_{i}-x_{j}}$ between all pairs of points $(x_{i},y_{i})$ and $(x_{j},y_{j})$. Under Assumptions ii), it is easy to see that

[TABLE]

which is closer to $\beta$ the larger $|\tilde{x}_{i}-\tilde{x}_{j}|$. In contrast, under Assumptions 2.2 we get that the expected slope between $y_{i,k}$ and $y_{j,l}$ equals (1) if $k\neq l$ and is zero for $k=l$. A similar argument holds for the median of slopes between and within groups. To include only the meaningful slopes for the estimation of the regression parameter $\beta$, we have to exclude all repeated measurement pairs within each group from the analysis.

Definition 2.6** (Block-Passing-Bablok regression for grouped data).**

For the estimation of the regression parameters $\alpha$ and $\beta$, compute the slopes of the connecting lines between any pair of points from different groups. The slopes are given by

[TABLE] 2. 2.

Discard identical measurements as well as all slopes with $S_{ij}^{kl}=-1$. Thereby obtain $N\leq\binom{n}{2}-\sum\limits_{k=1}^{m}\binom{p_{k}}{2}$ slopes. 3. 3.

The slope parameter $\beta$ is estimated by the shifted median of the slopes with the offset

[TABLE]

I.e. if the ranked sequence of slopes is given by $S_{(1)}\leq S_{(2)}\leq\dots\leq S_{(N)}$, $\beta$ is estimated by

[TABLE] 4. 4.

The intercept parameter $\alpha$ is estimated by

[TABLE]

Remark 2.7**.**

Naturally, the description of the groupwise method in Definition 2.6 is close to the description of the original method introduced by Passing and Bablok. If we set $m=n$ and $p_{k}=1$ we regain the classical estimator for the regression parameter $\beta$ [8]. In this case, we compute the slopes of the connecting lines between any pair of points, that are given by

[TABLE]

Without considering groups there are $\binom{n}{2}$ possible lines to connect any two points of an $n$-dimensional data set. Again two identical measurements with $x_{i}=x_{j}$ and $y_{i}=y_{j}$ are not considered for the estimation. Further, any slopes with a value of $-1$ are disregarded, such that we have at most $\binom{n}{2}$ slopes to consider.

Remark 2.8**.**

Note that setting $K=0$ in Definition 2.6 results in the Theil-Sen estimator. To get an estimator where methods can be used interchangeably, Passing and Bablok defined the offset determined by $K$. The definition of $K$ as the number of slopes with a value smaller than $-1$ in equation (2) corresponds to the null hypothesis $\beta=1$, and needs to be adapted for other hypothesis for the value of $\beta$. If the null hypothesis is not true, setting $K$ in this way introduces a bias towards higher estimates for $\beta$. The median slope between independent measurements within one group with offset $K$ is no longer given by zero if $K>0$. This effect can be seen in Table 1, e.g. for $\beta=0.2$ for both the classic and the Block-PBR. Since the offset $K$ drives the median of meaningless slopes between repeated measurements towards $1$, this can also lead to overconfidence for $\beta=1$ if the classic PBR is used. Table 1 for $\beta=0.8$ illustrates this effect.

3 Results

Definition 3.1**.**

Under Assumptions 2.2 let $\tilde{C}$ be the difference between the number of slopes that are larger than $\beta$ and the number of slopes lower than $\beta$. I.e. for

[TABLE]

[TABLE]

Note that we will denote this number by $C$ instead of $\tilde{C}$ in the special case of Assumptions ii) where no groups are considered.

Remark 3.2**.**

Using assumption 2.2 ii) we can write

[TABLE]

where $z_{i,j}^{k,l}=(\eta_{i,k}-\eta_{j,l})\beta^{-1}$ and $z_{i,j}^{{}^{\prime}k,l}=\epsilon_{i,k}-\epsilon_{j,l}$ are independent and from the same distribution [8]. Thus the number of slopes above and below $\beta$ and consequently $\tilde{C}$ are free from the regression parameters.

Theorem 1** (Variance of $\tilde{C}$).**

*If Assumptions 2.2 (i)-(ii) hold and $q_{ku}$ is the expected fraction of triplets with one point from group $u$ and two points from group $k$ with $\min(x_{i,k},x_{j,k})<x_{s,u}<\max(x_{i,k},x_{j,k})$ and

(a) if the group sizes are given by $(p_{1},\dots,p_{m})$,

[TABLE]

(b) Consequently, if the group sizes are equal, i.e. $p_{k}=p:=\tfrac{n}{m}$ for $k=1,\dots,m$,

[TABLE]

Remark 3.3** (classical Passing-Bablok regression).**

If we set $p=1$ and $m=n$ in Theorem 1 we regain the result from Passing and Bablok [8] and Theil-Sen [10] for independent measurements, where

[TABLE]

Remark 3.4**.**

In Theorem 1, note that (b) is a direct consequence from (a) as setting $p_{k}=p$ in equation (4) gives

[TABLE]

which leads to (5) as $mp=n$.

Theorem 2** (asymptotic normality of $\tilde{C}$).**

Let the number of groups $m$ be fixed and consider the limit of a large samplesize $n=\sum\limits_{i=1}^{m}p_{k}$.

(a) If the group sizes are given by $(p_{1}(n),\dots,p_{m}(n))$ and Assumptions 2.2 (i)-(ii) hold, $\tilde{C}$ is asymptotically normally distributed with mean zero and variance given by

[TABLE]

[TABLE]

where $q_{ku}$ is the expected fraction of triplets with one point from group $u$ and two points from group $k$ with $\min(x_{i,k},x_{j,k})<x_{s,u}<\max(x_{i,k},x_{j,k})$.

(b) More specifically, for $m$ groups with equal group sizes as in Theorem 1(b) and if Assumptions 2.2 and 2.3 hold, $\tilde{C}$ is asymptotically normally distributed with mean zero and variance

[TABLE]

Corollary 3.5** (confidence interval for $\beta$).**

Let $w_{\frac{\gamma}{2}}$ denote the \big{(}1-\frac{\gamma}{2}\big{)}-quantile of the standardized normal distribution. Let further $\tilde{\sigma}:=(\mathbb{V}[\tilde{C}])^{\frac{1}{2}}$,

[TABLE]

The asymptotic confidence interval for $\beta$ with significance level $1-\gamma$ is then given by

[TABLE]

Remark 3.6** (confidence interval for $\alpha$).**

Let $b_{L}$ denote the lower and $b_{U}$ denote the upper limit of the confidence interval $I$ for the slope $\beta$ in equation (6). Then,

[TABLE]

are the corresponding limits for the intercept $\alpha$.

Remark 3.7** (Equivalence of two methods with repeated measurements).**

Given two measurement methods with measurements given by $x_{i,k}$ and $y_{i,k}$.

• •

The equivalence of both methods can be concluded if $0\in[a_{L},a_{U}]$ and $1\in[b_{L},b_{U}]$.

• •

An intercept $\hat{\alpha}\neq 0$ indicates a constant systemic difference between the two measurement methods.

• •

A regression parameter $\hat{\beta}\neq 1$ indicates a proportional systemic difference between the measurement methods.

Note that here the asymptotic confidence level is calculated for $\beta$ and not $\alpha$.

Remark 3.8** (non-overlapping groups).**

In case of a data set where Assumption 2.3 is fulfilled, i.e. with non-overlapping groups on the x-axis, we can set $q_{ku}=0$ and formula (4) transforms to

[TABLE]

And consequently for equal group sizes

[TABLE]

Note that (9) equals the correction for tied ranks in one variable [10]. To see this consider the fact that setting $x_{i,k}$ to the mean of the corresponding group for $k=1,\dots,m$ does not change the sign of slopes between separated groups. So compared to the variance for non-overlapping groups the variance for overlapping groups is always smaller. For overlapping groups the test for the hypotheses $\beta=1$ from Remark 3.7 based on equation (9) is thus a conservative test for the equivalence of the measurement methods, but the power of the test could in principle be improved if a reliable estimate for $q_{ku}$ is available.

4 Discussion and empirical comparisons of

regular and Block-Passing-Bablok Regression

Passing-Bablok Regression is recommended as a robust method such that extreme values can be included and the errors do not have to be normally distributed [3]. However, the assumption of independent measurements in the classical Passing-Bablok regression seems to be frequently not fulfilled. Duplicated and repeated measurements are often an important part of studies to assess the variance of measurements and help to identify outliers [3]. As a random example among various studies including repeated measurements consider Figure 2 in [12], where measurements of different patients have been repeated in different numbers. But to which degree does this harm the results of a classical Passing-Bablok regression?

We recall, the Block-Passing-Bablok Regression only considers slopes between pairs of points from different groups for estimating the regression parameter $\beta$. Since the regression parameters are estimated as medians, a relatively small difference for the number of considered slopes does not have a large influence on the estimations. Hence if the group sizes are equal, instead of $\binom{n}{2}$ slopes as in the original method, we utilize approximately $\binom{n}{2}(1-\frac{1}{m})$ slopes. As shown in Theorem 2 (a) the discrepancy between the unadapted and the Block-Passing-Bablok Regression will quickly vanish for a large number of separated and equally sized groups.

In a hypothetical example with varying group sizes where $p_{1}\gg p_{2},\dots,p_{m}$, the Block-Passing-Bablok Regression yields better results than the original method. The median of all slopes between a very large and a very small group would be a slope within the large group. Since we assume random errors within groups, this median slope would be sampled from a set of meaningless slopes with mean [math], which transforms into $1$ if the offset $K$ is used. This biases the estimate of $\beta$ towards $1$ and lowers the power to detect $\beta\neq 1$. The Block-Passing-Bablok Regression provides more reliable estimates in this setting.

To illustrate the influence of the group sizes and the overlap between groups on the estimation of $\beta$ and the power of the test of the hypotheses $\beta=1$ in Remark 3.7 we simulated 16 illustrative scenarios. Evaluations of the Passing-Bablok Estimates have been performed with the mcr package [5] for the classic PBR and with a custom script based on mcr for the Block-PBR. The script is available from the authors upon request. A subset of the considered settings is illustrated in Figure 1. The results of the simulations are shown in Table 1.

5 Proofs

For Theorem 1, Theorem 2, and the verification of confidence intervals in Corollary 3.5 we need to compute the variance of

[TABLE]

the number of slopes $S_{ij}^{kl}$ bigger than $\beta$ minus the number of slopes $S_{ij}^{kl}$ smaller than $\beta$, as introduced in Definition 3.1.

To simplify the notation for the proofs, we will instead consider a transformed dataset, where $x_{i,k}$ are divided by $\beta$ such that the errors in both dimensions are identically distributed (see Assumption 2.2(ii)). If we also substract $\beta^{-1}x_{i,k}$ from $y_{i,k}$, the sign of the slopes in the transformed datasets determines whether a pair of points contributes to $P(\beta)$ or $Q(\beta)$ from Definition 3.1. In particular, all figures in the proof section are plotted with transformed coordinates, i.e. $(x_{i,k}^{\prime},y_{i,k}^{\prime})=(\beta^{-1}x_{i,k},y_{i,k}-\beta^{-1}x_{i,k})$.

5.1 Proof of Theorem 1

Proof.

By defining

[TABLE]

as the sign of the slope in the transformed dataset, it follows that

[TABLE]

We will hence compute

[TABLE]

Since $\mathbb{E}(Z_{ij}^{kl})=0$ for all \big{(}(i,k),(j,l)\big{)}, the second term on the right hand side vanishes and we have

[TABLE]

We split the right hand side into the sums over the squares $(Z_{ij}^{kl})^{2}$ and the sum over the remaining mixed terms $Z_{ij}^{kl}\cdot Z_{rs}^{tu}$:

[TABLE]

(a) The first term on the right hand side is not hard to compute since $(Z_{ij}^{kl})^{2}=1$ for all \big{(}(i,k),(j,l)\big{)} with $k\neq l$ holds. We get

[TABLE]

(b) For the computation of the mixed terms, we only need to consider those $Z_{ij}^{kl}\cdot Z_{rs}^{tu}$ with a common pair of indices, i.e. one of the four cases $(i,k)=(r,t),(i,k)=(s,u),(j,l)=(r,t)$ and $(j,l)=(s,u)$. Ignoring group assignments there are $4n(n-1)(n-2)$ such combinations. All other combinations vanish in expectation, due to independence. For the expectation of the sum of mixed terms, we look at the covariances of the $Z_{ij}^{kl}$ with each other. Since all cases lead to the same result, we can assume without loss of generality $(j,l)=(r,t)$ and multiply the result by four.We will combine triplets of mixed terms to simplify the calculations. By distinction of cases as sketched in Figure 2, we obtain for the expectation

[TABLE]

considering that the cases occur with probability $\frac{1}{3}$ and $\frac{2}{3}$, respectively.

Furthermore, we need to consider that we counted each pair of slopes three times due to the arrangement in triples. Consequently, we have to divide by $3$ and get the interim result

[TABLE]

So far we did not consider the groups. Since the method from Definition 2.6 does not consider any pairs of data points within the same group, we have to subtract those triples of indices \big{(}(i,k),(j,l),(s,u)\big{)} with $k=l,k=u$, and $l=u$, as well as $k=l=u$.

(c) There are

[TABLE]

possible combinations of three elements from two different groups. Let’s say $k=l$. In that case, the only remaining term in (13) is $Cov(Z_{is}^{ku},Z_{js}^{lu})$ and we have to subtract the other two terms from our computation of (11) because we wrongfully included them in step (b). In case of strictly separated groups, however, those combinations vanish as can be seen below in equation (16).

In case of two equal indices there are three different possibilities, how the points are located relativ to each other. Looking at the single point from group $u$, it can be located above, in between or below the two points from group $k$ (relating to the y-axis). Those arrangements occur with probability $\frac{1}{3}$ each and are illustrated in Figure 3 if $x_{s,u}>\max(x_{i,k},x_{j,k})$ and in Figure 4 for $\min(x_{i,k},x_{j,k})<x_{s,u}<\max(x_{i,k},x_{j,k})$. Summed up, those arrangements yield an expectation of

[TABLE]

(d) Further, in (14) we wrongfully included the case $k=l=u$. This means, we have to substract

[TABLE]

and can combine (12) and (14)-(17) to get for the term in (11)

[TABLE]

the wanted formula for the variance. Naturally, this formula holds likewise for the less complex settings of equally sized non-overlapping groups with $p_{k}=p$ for all $k$, i.e.

[TABLE]

and for non-grouped data with $m=n$ and $p_{k}=1$ we regain the classic result

[TABLE]

5.2 Proof of Theorem 2

Proof.

The proofs for the PBR and TSR for non-grouped measurements [8, 10] relied on results for rank correlation methods. Our proof will follow the same route as shown in chapter 5 of Kendalls book [4]. To show the asymptotic normality we will show that the moments of the distribution of $\tilde{C}$ tend to those of the normal distribution [1, Thm 30.2].

Under the null hypothesis we have that $\mathbb{P}(Z^{i,k}_{j,l}=1)=\mathbb{P}(Z^{i,k}_{j,l}=-1)$. Hence moments of $\tilde{C}$ of odd order vanish. For moments of even order, we need to compute

[TABLE]

Consider the expansion of the sum above, which consists of summands with $2r$ factors each. Since $\mathbb{E}[Z_{ij}^{kl}]=0$ and the $Z_{ij}^{kl}$ are independent if $k\neq l$ each summand with an independent factor will vanish in the expectation. Let us consider the number of summands where the factors are pairwise linked by exactly one suffix. Each of these summands will look like this:

[TABLE]

which due to independence could be split up into the form

[TABLE]

However, by doing so we would no longer connect the cases as in equation (13). So we ignore the above and combine cases, as we have done in the proof of Theorem 1:

For each $Z_{ij}^{kl}Z_{is}^{ku}$ in each summand as shown in equation (22) there are two other combinations that are exactly the same except that there is an $Z_{ij}^{kl}Z_{jl}^{us}$ or $Z_{is}^{ku}Z_{jl}^{us}$ instead of $Z_{ij}^{kl}Z_{is}^{ku}$, respectively. Therefore we define for each triplet of points $t=((i,k),(j,l),(s,u))$ the pairings

[TABLE]

For each summand with pairwise tied indices let $t_{1},\dots,t_{r}$ be the triplet of indices used in the $r$-th pairing. Let us then consider the set of $3^{r}$ sums which is represented by $\prod\limits_{i=1}^{r}\sum\limits_{j=1}^{3}p_{j}(t_{i}).$ Now we see that

[TABLE]

if the points of $t_{1},\dots,t_{r}$ are a disjoint and hence independent set of $3r$ points.

The remaining questions are:

1. (a)

How many summands with pairwise tied indices are there in the expanded sum ignoring group associations? 2. (b)

Which of these need to be omitted due to group associations for non-overlapping groups? 3. (c)

What changes for overlapping groups?

Answer to (a):

There are $\binom{2r}{r}$ ways to choose the first factors of $r$ pairings and $r!$ ways to assign the remaining $r$ factors. But now we have counted some combinations twice so we have to divide by $2^{r}$ to get the final number of ways to choose $r$ pairings. Given such a combination there are now $3r$ different indices to choose, so $\frac{n!}{(n-3r)!}\sim n^{3r}$ possibilities.

Thus we end up with

[TABLE]

combinations with $r$ pairs with tied indices. Finally since we counted each summand $3^{r}$ times we would get for (21)

[TABLE]

if we ignore group associations.

Answer to (b):

For $t=((i,k),(j,l),(s,u))$ there are 3 different scenarios.

Scenario 1: all groups are different ($k\neq l\neq u$)

No difference to the above calculation

Scenario 2: all groups are equal ($k=l=u$)

The corresponding summand is not included in (21) which we mimic by setting

\sum\limits_{j=1}^{3}\mathbb{E}\Big{[}p_{j}(t)\Big{]}=0

Scenario 3: two groups are equal ($k=l\neq u$ or $k\neq l=u$ or $k=u\neq l$)

Althoug some summands are not included in this case, we do not have to change anything since the effect vanishes for separated groups as shown in equation (16).

For separated groups only scenario 2 alters the above calculation. So we have to subtract the number of combinations with at least one factor where all indices are from the same group. The probablity to choose 3 indices from group $k$ is given by $n^{-3}p_{k}(n)(p_{k}(n)-1)(p_{k}(n)-2)$ and thus the fraction of all sums that have no $t_{i}$ from scenario 2 tends to

[TABLE]

which simplifies for $m$ groups with equal sizes to

[TABLE]

So we have to multiply equation (23) and (24) and end up with

[TABLE]

Answer to (c):

If groups can overlap, let $r_{1}$ be the number of triplets where groups do not overlap (Scenario 3.1) and $r_{2}$ the number of triplets where groups overlap (Scenario 3.2).

In this case we get

[TABLE]

So for each summand with an odd number of triplets with $\min(x_{i,k},x_{j,k})<x_{s,u}<\max(x_{i,k},x_{j,k})$ we have to substract $\frac{2}{3^{r}}$ just as in the proof of Theorem 1. We set

[TABLE]

The probability to get a summand with odd $r_{2}$ given there is no triplet from scenario 2 is now given by $\frac{1}{2}(1-(1-2\frac{b}{a})^{r}).$ This leads us to

[TABLE]

As all other terms are of lower order in $n$ this completes the proof. To see this consider the number of summands where the factors are all linked but not pairwise, i.e. there is at least one summand that is linked to more than one other summand. In this case, there are less than $3r-1$ indices to choose, such that the order is not greater than $n^{3r-1}$.

5.3 Justification of confidence intervals in Corollary 3.5

For the justification of the formula $M_{1}=\lfloor\frac{N-\tilde{C}_{\gamma}}{2}\rfloor$ we proceed just as Passing and Bablok [8]. Consider

[TABLE]

which holds if and only if $M_{1}+K\leq~{}Q(\beta)$ and $M_{1}-K\leq~{}P(\beta)$ hold and thus if and only if

[TABLE]

We know that the distribution of $\tilde{C}:=P(\beta)-Q(\beta)$ does not depend on the distribution of $(X,Y)$, whereas the distribution of $K$ does. Therefore, we can not provide a general formula for $M_{1}$ satisfying

[TABLE]

that does not dependent on the distribution of $(X,Y)$. However, as shown in Theorem 2, $\tilde{C}$ is asymptotically normally distributed, such that P\big{(}-\tilde{C}_{\gamma}\leq\tilde{C}\leq\tilde{C}_{\gamma}\big{)}\sim 1-\gamma. Thus, we argue just as Passing and Bablok [8] that $M_{1}$ can be found by setting $N-2M_{1}=\tilde{C}_{\gamma}$.

acknowlegement

We thank Peter Pfaffelhuber for various discussions pointing us in the right direction as well as Lukas Steinberger and an anonymous reviewer for helpful comments and careful reading of the manuscript.