Die Rohdaten der Studie mit dem Titel “Test-retest reliability and minimal detectable change scores for the short physical performance battery, one-legged standing test and timed up and go test in patients undergoing hemodialysis” finden Sie auf Moodle. Auf den Volltext können Sie hier zugreifen.

Übung 1

Aufgaben

  1. Laden Sie die Daten in R. Splitten Sie den Datensatz auf nach Test. Erstellen Sie dazu jeweils neue R-Objekte. Bringen Sie dann die drei Teildatensätze in ein long-format. Sie können für den zweiten Schritt folgenden Code verwenden:
library(tidyverse)
df.long <- df %>% 
  pivot_longer(cols = c(X,Y), values_to = "value_name", names_to = "time_point")

  1. Reproduzieren Sie die angegebenen SEM in Tabelle 4. Verwenden Sie das Paket SimplyAgree (siehe Slides für den R-Code).

  2. Ohne auf die Details des ICCs (Tabelle 3) einzugehen: Haben die Autor:innen den SEM mit oder ohne Berücksichtigung eines systematischen Fehlers berechnet?

Lösungen

Hier wird der Datensatz importiert, gesplitted nach Test und dann jeweils ins long-format gebracht. Natürlich gibt es auch andere Möglichkeiten, dies zu tun.

library(rio)
df.tests <- import("../data/data-villar.csv")

library(tidyverse)
df.sppb <- df.tests %>% 
  select(ID, SPPB_A, SPPB_B)

df.sppb.long <- df.sppb %>% 
  pivot_longer(cols = c(SPPB_A, SPPB_B), values_to = "value", names_to = "time")

df.olst <- df.tests %>% 
  select(ID, OLST_A, OLST_B)

df.olst.long <- df.olst %>% 
  pivot_longer(cols = c(OLST_A, OLST_B), values_to = "value", names_to = "time")

df.tug <- df.tests %>% 
  select(ID, TUG_A, TUG_B)

df.tug.long <- df.tug %>% 
  pivot_longer(cols = c(TUG_A, TUG_B), values_to = "value", names_to = "time")

Die drei SEM erhalten wir wie folgt:

SPPB:

library(SimplyAgree)

reli_stats(
  measure = "value",   # name of the column with ratings
  item = "time",      # rater or time identifier
  id = "ID",      # subject identifier
  data = df.sppb.long, # your data frame
  se_type = "ICC2",   # this is important!
  cv_calc = "MSE",     # how to calculate CV
  conf.level = 0.95,
  other_ci = TRUE
)
## 
## Coefficient of Variation (%): 7.04 95% C.I. [6.04, 8.45]
## Standard Error of Measurement (SEM): 0.71 95% C.I. [0.61, 0.851]
## Standard Error of the Estimate (SEE): 0.69 95% C.I. [0.592, 0.827]
## Standard Error of Prediction (SEP): 0.99 95% C.I. [0.85, 1.19]
## 
## Intraclass Correlation Coefficients
##            Model         Measures  Type    ICC Lower CI Upper CI
## 1 one-way random        Agreement  ICC1 0.9442   0.9183   0.9620
## 2 two-way random        Agreement  ICC2 0.9442   0.9178   0.9622
## 3  two-way fixed      Consistency  ICC3 0.9459   0.9206   0.9632
## 4 one-way random   Avg. Agreement ICC1k 0.9713   0.9574   0.9806
## 5 two-way random   Avg. Agreement ICC2k 0.9713   0.9571   0.9808
## 6  two-way fixed Avg. Consistency ICC3k 0.9722   0.9587   0.9813

Mit “ICC2” erhalten Sie den SEM mit Berücksichtigung des systematischen Fehlers, mit “ICC3” ohne. Bis auf die zweite Stelle nach dem Komma stimmt das Resultat.

Die manuelle Berechnung mit der klassischen Formel (wie in der Methodik angegeben) via SD und ICC führt fast zu einem identischen Resultat:

sd(c(df.sppb$SPPB_A, df.sppb$SPPB_B), na.rm = TRUE) * sqrt(1-0.9442)
## [1] 0.6904926

Warum die Autor:innen auf 0.72 kommen, weiss ich auch nicht. Es wird dort nicht beschrieben, welche SD sie nehmen. Hier wird die Standardabweichung über beide Messungen genommen.

Auf Grund der Beschreibung in der Methode ist anzunhemen, dass die Autor:innen den systematischen Fehler berücksichtigt haben.

“two-way random-effects model were used to assess the test–retest reliability of the data”


OLST:

library(SimplyAgree)

reli_stats(
  measure = "value",   # name of the column with ratings
  item = "time",      # rater or time identifier
  id = "ID",      # subject identifier
  data = df.olst.long, # your data frame
  se_type = "ICC2",   # this is important!
  cv_calc = "MSE",     # how to calculate CV
  conf.level = 0.95,
  other_ci = TRUE
)
## 
## Coefficient of Variation (%): 34.3 95% C.I. [28.9, 42.4]
## Standard Error of Measurement (SEM): 4.99 95% C.I. [4.27, 6.01]
## Standard Error of the Estimate (SEE): 4.71 95% C.I. [4.03, 5.67]
## Standard Error of Prediction (SEP): 6.86 95% C.I. [5.87, 8.26]
## 
## Intraclass Correlation Coefficients
##            Model         Measures  Type    ICC Lower CI Upper CI
## 1 one-way random        Agreement  ICC1 0.8905   0.8405   0.9255
## 2 two-way random        Agreement  ICC2 0.8905   0.8405   0.9255
## 3  two-way fixed      Consistency  ICC3 0.8905   0.8403   0.9256
## 4 one-way random   Avg. Agreement ICC1k 0.9421   0.9134   0.9613
## 5 two-way random   Avg. Agreement ICC2k 0.9421   0.9134   0.9613
## 6  two-way fixed Avg. Consistency ICC3k 0.9421   0.9132   0.9614

Hier sind ICC2 und ICC3 gleich (kein Bias). Die Abweichung zum publizierten Wert ist aber etwas grösser.

Auch mit der manuellen Berechnugn sind die Werte nicht exakt gleich.

sd(c(df.olst$OLST_A, df.olst$OLST_B), na.rm = TRUE) * sqrt(1-0.8905)
## [1] 4.928476


TUG:

library(SimplyAgree)

reli_stats(
  measure = "value",   # name of the column with ratings
  item = "time",      # rater or time identifier
  id = "ID",      # subject identifier
  data = df.tug.long, # your data frame
  se_type = "ICC2",   # this is important!
  cv_calc = "MSE",     # how to calculate CV
  conf.level = 0.95,
  other_ci = TRUE
)
## 
## Coefficient of Variation (%): 10.6 95% C.I. [9.11, 12.8]
## Standard Error of Measurement (SEM): 1.18 95% C.I. [1.01, 1.41]
## Standard Error of the Estimate (SEE): 1.16 95% C.I. [0.995, 1.39]
## Standard Error of Prediction (SEP): 1.65 95% C.I. [1.42, 1.98]
## 
## Intraclass Correlation Coefficients
##            Model         Measures  Type    ICC Lower CI Upper CI
## 1 one-way random        Agreement  ICC1 0.9644   0.9476   0.9758
## 2 two-way random        Agreement  ICC2 0.9644   0.9476   0.9759
## 3  two-way fixed      Consistency  ICC3 0.9645   0.9477   0.9760
## 4 one-way random   Avg. Agreement ICC1k 0.9819   0.9731   0.9878
## 5 two-way random   Avg. Agreement ICC2k 0.9819   0.9731   0.9878
## 6  two-way fixed Avg. Consistency ICC3k 0.9819   0.9731   0.9878

Auch hier sind ICC2 und ICC3 (praktisch) gleich (kein Bias). Die Abweichung zum publizierten Wert ist aber etwas grösser.

Auch mit der manuellen Berechnugn sind die Werte nicht exakt gleich.

sd(c(df.tug$TUG_A, df.tug$TUG_B), na.rm = TRUE) * sqrt(1-0.9644)
## [1] 1.14371

Das 95% CI in der Publikation sieht etwas seltsam aus…

Übung 2

Aufgaben

  1. Reproduzieren Sie auch die angegebenen SDCs.
  2. Rechnen Sie zusätzlich einen \(SDC_{95}\)
  3. Interpretieren Sie die Werte.

Lösungen

Wir bestimmen zuerst die z-Werte:

z.9 <- qnorm(0.95)
z.95 <- qnorm(0.975)

Damit können wir nun die SDCs bestimmen. Hier werden die SEM-Werte aus der Studie verwendet, man könnte aber auch jene von oben nehmen:

(SDC.9.sppb <- z.9 * sqrt(2) * 0.72)
## [1] 1.674846
(SDC.95.sppb <- z.95 * sqrt(2) * 0.72)
## [1] 1.995702
(SDC.9.olst <- z.9 * sqrt(2) * 4.82)
## [1] 11.21216
(SDC.95.olst <- z.95 * sqrt(2) * 4.82)
## [1] 13.36011
(SDC.9.tug <- z.9 * sqrt(2) * 1.24)
## [1] 2.884456
(SDC.95.tug <- z.95 * sqrt(2) * 1.24)
## [1] 3.437041

Wenn wir mit dem entsprechenden Assessment bei einer Person der entsprechenden Population eine Differenz von z.B. 1.7 (SPPB) messen, können wir zu 90% sicher sein, dass es eine wahre Veränderung und nicht nur Messfehelr ist. (bzw. 95% sicher sein, wenn die Differenz \(\ge 2\) ist. Analog für die anderen Assessments.

Übung 3

Aufgaben

  1. Untersuchen Sie für das Assessment SPPB den systematischen Fehler. Berechnen Sie dafür zuerst die durchschnittliche Differenz zwischen den zwei Messzeitpunkten.

  2. Erstellen Sie einen Bland-Altman Plot mit mit einem Agreement-Level von 95% (Syntax, siehe Slides). Vergleichen Sie mit Abbildung 2 aus der Publikation.

  3. Interpretieren Sie den Bias und die LoA.

Lösungen

Die durchschnittliche Differenz berechnet sich wie folgt:

mean(df.sppb$SPPB_A - df.sppb$SPPB_B, na.rm = TRUE)
## [1] -0.2

Für die Bland-Altman-Analyse verwenden wir das Paket “SimplyAgree”:

ba.sppb <- agree_test(x = df.sppb$SPPB_A,
                 y = df.sppb$SPPB_B,
                agree.level = .95)

plot(ba.sppb)

Der Plot sieht genau gleich aus (schöner 🤓) wie jener in der Studie.

Der Analyse entnehmen wir den oben schon manuell berechneten Bias. Da das 95% CI die 0 knapp beinhaltet, unterscheidet sich dieser nicht signifikant von 0.

print(ba.sppb)
## Limit of Agreement = 95%
## 
## ###- Shieh Results -###
## Exact 90% C.I.  [-2.3423, 1.9423]
## Hypothesis Test: No Hypothesis Test
## 
## ###- Bland-Altman Limits of Agreement (LoA) -###
##           Estimate Lower CI Upper CI CI Level
## Bias        -0.200  -0.4407  0.04072     0.95
## Lower LoA   -2.104  -2.4495 -1.75858     0.90
## Upper LoA    1.704   1.3586  2.04950     0.90
## 
## ###- Concordance Correlation Coefficient (CCC) -###
## CCC: 0.9416, 95% C.I. [0.9063, 0.9638]

Das könnte man auch mit einem einfachen t-Test erfahren:

t.test(df.sppb$SPPB_A - df.sppb$SPPB_B)
## 
##  One Sample t-test
## 
## data:  df.sppb$SPPB_A - df.sppb$SPPB_B
## t = -1.6598, df = 64, p-value = 0.1018
## alternative hypothesis: true mean is not equal to 0
## 95 percent confidence interval:
##  -0.44071797  0.04071797
## sample estimates:
## mean of x 
##      -0.2

Warum in der Studie dafür ein P-Wert von 0.942 angegeben ist (siehe hier), konnte ich bis jetzt nicht herausfinden.