3 IVD检测方法一致性对比

3.1 模拟数据

```{sas}
data voriconazole;
  call streaminit(12345); 
  do id = 1 to 130; 
    conc_x = rand("Normal", 2.0, 0.9);
      conc_y = conc_x + rand("Normal", -0.1, 0.3);
    output;
  end;
  drop id;
run;

proc print data=voriconazole(obs=10);run;

proc sgplot data=voriconazole;
  scatter x=conc_x y=conc_y;
run;
```

首先模拟生成130个样本，每个样本包含血浆药浓度检测值conc_x和conc_y。conc_x为LC-MS/MS检测结果，conc_y为清谱科技（苏州）有限公司生产的伏立康唑药物浓度检测试剂盒（原位电离质谱法）及便携式质谱分析系统检测结果。conc_x的均值和标准差分别设定为2.0 $\mu g/mL$和0.9$\mu g/mL$。conc_y的均值和标准差分别设定为1.9 $\mu g/mL$和0.95$\mu g/mL$. 使用以上SAS代码输出前10组样本和样本散点图如下：

3.2 散点图

3.2.1 Bland-Altman 一致性界限的置信区间

一致性界限的置信区间为 \[ \text{Lower LOA CI} = \text{Lower LOA} \pm z_{1-\alpha/2}\cdot SE_{LOA} \] \[ \text{Upper LOA CI} = \text{Upper LOA} \pm z_{1-\alpha/2}\cdot SE_{LOA} \]

如果样本量较少t distribution 不能近似成正态分布。需要将上文中的$z_{1-\alpha/2}$ 或$z_{1-(1-agree)/2}$ 替换成$t_{1-\alpha/2,\ df}$或$t_{1-(1-agree)/2, df}$

3.2.2 绝对偏倚Bland-Altman

先用data步计算每个样本间conc_x和conc_y的均值和差，保存为intra_mean和intra_diff。接下来用proc means计算intra_diff（每个样本中考核试剂与 LC-MS/MS 检测血浆测定值差值）的均值和标准差。再用proc sql计算95%一致性界限并将其保存为lower和upper。最后使用proc sgplot绘制散点图，用refline绘制水平线。

```{sas}
data voriconazole_diff;
    set voriconazole;
    intra_mean = (conc_y + conc_x)/2;
    intra_diff = conc_y - conc_x;
run;
proc means data=voriconazole_diff noprint;
  var intra_diff;
  output out=limits mean=mean_diff std=sd_diff;
run;
proc sql noprint ;
  select mean_diff-1.96*sd_diff, mean_diff+1.96*sd_diff, mean_diff
  into   :lower,  :upper ,:mean_diff
  from limits ;
quit;

proc sgplot data=voriconazole_diff;
  scatter x=intra_mean y=intra_diff / markerattrs=(symbol=CircleFilled color=black);
  refline 0 /lineattrs=(pattern=shortdash color=gray);
  refline &mean_diff &lower &upper;
  xaxis label="测定均值" labelattrs=(family="SimSun" size=12);
  yaxis label="测定值之差" labelattrs=(family="SimSun" size=12);
run;
```

结果如下：

图中三条实心横线从上到下分别为：考核试剂与 LC-MS/MS 检测血浆测定值差值的95%一致性界限上限；考核试剂与 LC-MS/MS 检测血浆测定值差值的均值；考核试剂与 LC-MS/MS 检测血浆测定值差值的95%一致性界限下限。虚线是考核试剂与 LC-MS/MS 检测血浆测定值差值为0。

3.2.3 相对偏倚Bland-Altman

与上一节类似，只不过这次我们计算每个样本间conc_x和conc_y的均值和比值。

```{sas}
data voriconazole_ratio;
    set voriconazole;
    intra_mean = (conc_y + conc_x)/2;
    intra_ratio = conc_y /conc_x;
run;
proc means data=voriconazole_ratio noprint;
  var intra_ratio;
  output out=limits mean=mean_ratio std=sd_ratio;
run;
proc sql noprint ;
  select mean_ratio-1.96*sd_ratio, mean_ratio+1.96*sd_ratio, mean_ratio
  into   :lower,  :upper ,:mean_ratio
  from limits ;
quit;

proc sgplot data=voriconazole_ratio;
  scatter x=intra_mean y=intra_ratio / markerattrs=(symbol=CircleFilled color=black);
  refline 0 /lineattrs=(pattern=shortdash color=gray);
  refline &mean_ratio &lower &upper;
  xaxis label="测定均值" labelattrs=(family="SimSun" size=12);
  yaxis label="测定值之比" labelattrs=(family="SimSun" size=12);
run;
```

结果如下：

图中三条实心横线从上到下分别为：考核试剂与 LC-MS/MS 检测血浆测定值比值的95%一致性界限上限；考核试剂与 LC-MS/MS 检测血浆测定值比值的均值；考核试剂与 LC-MS/MS 检测血浆测定值比值的95%一致性界限下限。虚线表示比值为1。

3.2.4 直线回归散点图

使用proc reg进行简单线性回归，其中conc_x为自变量；conc_y为因变量。使用proc sgplot绘制线性回归散点图。

```{sas}
ODS OUTPUT ParameterEstimates=reg_params FitStatistics=reg_fit;
proc reg data=voriconazole;
  model conc_y = conc_x;
run;
ODS OUTPUT CLOSE;

proc sgplot data=voriconazole;
  reg x=conc_x y=conc_y / lineattrs=(color=blue thickness=2);
  scatter x=conc_x y=conc_y / markerattrs=(color=black symbol=CircleFilled);
  xaxis label="LC-MS/MS检测血浆定值";
  yaxis label="考核试剂测定血浆值";
run;
```

结果如下：

图中红框1为截距a，红框2为斜率b。截距和斜率分别为-0.06690和0.98406。回归模型散点图如下：

3.3 PROC CORR 计算相关性系数

3.3.1 Pearson相关系数

使用PROC CORR Pearson计算conc_x和conc_y相关性系数，var后面放需要计算相关性的变量。

```{sas}
proc corr data=voriconazole pearson;
  var conc_x conc_y;
run;

```

结果如下：

图中红框1为pearson相关性系数值0.93827。红框2为相关性系数值是否为0的检验，此处我们拒绝相关性系数值为0，接受相关性系数值不为0，为0.93827。

3.3.2 Spearman相关系数

使用PROC CORR spearman计算相关性系数，var后面放需要计算相关性的变量。

```{sas}
proc corr data=voriconazole spearman;
  var conc_x conc_y;
run;

```

结果如下：

得出spearman相关性系数值为0.93346。

3.4 医学决定水平附近的检测结果分析

根据公式： \[ \hat{B}_{c} = \hat{Y}_{c}-X_c \] \[ \hat{B}_{c}\text{ CI} = \hat{Y}_{c} \text{ CI}-X_c \]

推导：$Equation 3.2

编写代码如下

```{sas}
data newdata1;
    input conc_x;
    datalines;
    0.5
    5
;
run;
data newdata2;
 set voriconazole newdata1;
 run;

proc reg data=newdata2;
  id conc_x;
  model conc_y = conc_x/r clb clm;
  output out=ypredclm
    p=yhat zhat
    r=yresid zresid
    stdp=SE
    lclm=lclm
    uclm=uclm;
run;

```

```{sas}

data ypredclm2;
    set  ypredclm;
    where conc_x=0.5 or conc_x = 5;
run;

data bias;
    set ypredclm2;
    bias = yhat-conc_x;
    LCI = lclm-conc_x;
    UCI = uclm-conc_x;
    keep conc_x bias SE LCI UCI;
run;
```

结果如下：

以医学决定水平Xc_value为0.5和5为例， bias为预期偏倚，LCI和UCI分别为预期偏倚的置信区间下限和上限

3.5 使用R验证结果准确性

因为sas没有将这些计算公式封装成函数，所以我们在R中使用blandr包和mcr包来验证SAS结果的准确性。

library(blandr)
library(tidyverse)
library(openxlsx)
library(flextable)
library(mcr)

先把之前在sas里生成的模拟数据导入到R中，数据和sas保持一致。

df_simulation <- read.xlsx('Data/simulation.xlsx')
head(df_simulation)%>%flextable()%>%autofit()

conc_x	conc_y
2.237810	2.460228
2.736132	2.470299
3.386130	2.915984
1.872618	2.085219
2.059159	2.326737
1.866526	1.853563

3.5.1 绝对偏倚Bland-Altman

先看下在sas里面我们算出来的conc_x和conc_y的差值均值与95%置信区间。

我们使用blandr.statistics这个函数计算一致性界限上下限，结果和与SAS完全一致。

bland_diff <- blandr.statistics(df_simulation$conc_y,df_simulation$conc_x, sig.level=0.95)
bland_diff$bias

[1] -0.09802023

bland_diff$lowerLOA

[1] -0.6859971

bland_diff$upperLOA

[1] 0.4899566

3.5.2 线性回归

我们计算医学决定水平的预期偏倚及其95%置信区间也就是$\hat{B}_c, \hat{B}_{c,low}, \hat{B}_{c,high}$出自于 EP09A3-Measurement Procedure Comparison and Bias Estimation Using Patient Samples Approved Guideline-Third Edition EP09 Measurement Procedure Comparison and Bias Estimation Using Patient Samples (n.d.) 。我们接下来会用到的R mcr就是以此书作为主要参考文献。在SAS中我们使用proc reg进行线性回归，得到的斜率和截距分别为0.98406和-0.06690。我们在R中使用mcr包mcreg函数进行线性回归，可以看到斜率和截距与SAS结果一致。

lm_model <- mcreg(df_simulation$conc_x, df_simulation$conc_y, 
                  method.reg = "LinReg", method.ci = "analytical", 
                  mref.name = 'LC-MS/MS plasma', mtest.name = 'IVD plasma')
printSummary(lm_model)



------------------------------------------

Reference method: LC-MS/MS plasma
Test method:     IVD plasma
Number of data points: 130

------------------------------------------

The confidence intervals are calculated with analytical method.
Confidence level: 95%


------------------------------------------

LINEAR REGRESSION FIT:

                  EST         SE        LCI        UCI
Intercept -0.06690167 0.06791444 -0.2012820 0.06747865
Slope      0.98405515 0.03206456  0.9206099 1.04750036

NULL

做线性回归图

plot(lm_model)

3.5.3 医学决定水平附近的检测结果分析

我们使用MCResult.calcBias函数计算医学决定水平附近的预期偏倚及95%置信区间。这里我们以医学决定水平0.5，2.5，5为例，与sas结果差距很小。

calcBias(lm_model, x.levels = c(0.5, 5))

   Level       Bias         SE        LCI        UCI
X1   0.5 -0.0748741 0.05350567 -0.1807442 0.03099601
X2   5.0 -0.1466259 0.10124370 -0.3469539 0.05370203

3.6 Appendix

1 两方法检测结果之差的一致性界限及一致性界限的置信区间出自Bland 和Altman的论文 Bland and Altman (1986)。 \[ LOA = \bar{d}\pm z_{1-({1-agree}/{2})} \cdot S_d\] 其中$\bar{d}$为两方法检测结果差值平均值，$S_d$为两方法检测结果差值的标准差。首先我们是用Delta method从$Var(S_d^2)$推导出$Var(S_d)$。

让$\theta = S_d^2$

\[ \begin{align} Var(S_d) &= Var(\theta^{\frac{1}{2}}) \\ &= [g'(\theta)]^2\cdot Var(\theta) \\ &= [\frac{d\ \theta^{\frac{1}{2}}}{d\ \theta}]^2\cdot Var(\theta) \\ &= [\frac{1}{2}S_d^{-1}]^2\cdot \frac{2S_d^4}{n-1} \\ &= \frac{S_d^2}{2(n-1)} \end{align} \tag{3.1}\] LOA的variance可以计算为 \[ \begin{align} Var(LOA) &= Var(\bar{d}+ z \cdot S_d)\\ &= Var(\bar{d})+z^2Var(S_d)\\ &= S_d^2/n+z^2\cdot \frac{S_d^2}{2(n-1)} \\ &= S_d^2\left( \frac{1}{n} +\frac{z^2}{2(n-1)}\right) \end{align} \]

因此 \[SE_{LoA} = S_d \cdot \sqrt{\frac{1}{N}+ \frac{z^2_{1-(1-agree)/2}}{2 \cdot(N-1)} }\]

\[ \begin{aligned} \hat{B}_{c} \text{ CI} &= \hat{B}_c \pm 2 S_{y \cdot x} \sqrt{\frac{1}{N} + \frac{(X_c - \bar{X})^2}{\sum(X_i - \bar{X})^2}}\\ &=\hat{Y}_{c}-X_c\pm 2 S_{y \cdot x} \sqrt{\frac{1}{N} + \frac{(X_c - \bar{X})^2}{\sum(X_i - \bar{X})^2}}\\ &=\hat{Y}_{c} \text{ CI}-X_c \end{aligned} \tag{3.2}\]

Reference

Bland, J. M., and D. G. Altman. 1986. “Statistical Methods for Assessing Agreement Between Two Methods of Clinical Measurement.” Lancet 1 (8476): 307–10.

EP09 Measurement Procedure Comparison and Bias Estimation Using Patient Samples. n.d. Accessed September 26, 2025. https://clsi.org/shop/standards/ep09/.