**```{r child = ‘descriptive.Rmd’, eval=TRUE}

**```

WP1

Investigating the effects of intervention on treatment group and comparing it with control group (double blind study)

# here is where the data was converted to long format, including it to confirm if it was done correctly

long_data_wmc <- wp1_data %>%
  pivot_longer(cols = c(bl_wmc_sum, fup6_wmc_sum, fup12_wmc_sum), names_to = "time_point", values_to = "outcome")

long_data_hhs <- wp1_data %>%
  pivot_longer(cols = c(bl_hhs_total, fup_hhs_total), names_to = "time_point", values_to = "outcome")

Analysis

Normality testing

Shapiro-Wilk normality testing

This is a summary of all the shapiro wilk tests for normality carried out below

baseline WOMAC in group 1 : normally distributed
baseline WOMAC in group 2 : normally distributed
baseline HHS in group 1 : normally distributed
baseline HHS in group 2 : not normally distributed
follow-up WOMAC in group 1 : not normally distributed
follow-up WOMAC in group 2 : not normally distributed
follow-up HHS in group 1 : not normally distributed
follow-up HHS in group 2 : not normally distributed

I used parametric methods for testing (based on large sample size), we need to confirm from a statistician whether this is okay considering the distribution of data.

## [1] "Shapiro-Wilk test for testing normality of distribution of baseline WOMAC in both groups "
## $`1`
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.97856, p-value = 0.5065
## 
## 
## $`2`
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.97911, p-value = 0.4889
## [1] "Histogram of distribution of baseline WOMAC in both groups "

## [1] "Shapiro-Wilk test for testing normality of distribution of baseline HHS in both groups "
## $`1`
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.96669, p-value = 0.1114
## 
## 
## $`2`
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.94825, p-value = 0.01777
## [1] "Histogram of distribution of baseline HHS in both groups "

## [1] "Shapiro-Wilk test for testing normality of distribution of follow-up WOMAC in both groups "
## $`1`
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.93397, p-value = 0.04547
## 
## 
## $`2`
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.84595, p-value = 2.885e-05
## [1] "Histogram of distribution of follow-up WOMAC in both groups "

## [1] "Shapiro-Wilk test for testing normality of distribution of follow- up HHS in both groups "
## $`1`
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.89649, p-value = 0.002032
## 
## 
## $`2`
## 
##  Shapiro-Wilk normality test
## 
## data:  X[[i]]
## W = 0.68373, p-value = 1.457e-08
## [1] "Histogram of distribution of follow-up HHS in both groups "

Results

T-tests

Comparing mean WOMAC in both groups at baseline and follow-up

Baseline: The t-test was conducted to compare the baseline WOMAC scores of Group 1 and Group 2. The calculated t-value was 0.328 with 99 degrees of freedom. The associated p-value was 0.7433, which is more than the chosen significance level of 0.05. Since the p-value is more than the chosen significance level (α = 0.05), we accept the null hypothesis. This indicates that there is no significant difference in baseline WOMAC scores between Group 1 and Group 2.

Follow-up: The t-test was conducted to compare the follow-up womac scores of Group 1 and Group 2. The calculated t-value was 7.813 with 76 degrees of freedom. The associated p-value was <0.001 , which is less than the chosen significance level of 0.05. Since the p-value is less than the chosen significance level (α = 0.05), we reject the null hypothesis. This indicates that there Was a significant difference in follow-up WOMAC scores between Group 1 and Group 2. The effect size for the difference was large (Cohen’s d > 0.8).

## WOMAC scores in treatment and control gorup at baseline
## 
##  Two Sample t-test
## 
## data:  bl_wmc_sum by tvc
## t = 0.32835, df = 99, p-value = 0.7433
## alternative hypothesis: true difference in means between group 1 and group 2 is not equal to 0
## 95 percent confidence interval:
##  -4.662974  6.512268
## sample estimates:
## mean in group 1 mean in group 2 
##        40.69388        39.76923
## 
##  Effect size (cohen's D)
## [1] 0.06537309
## WOMAC scores in treatment and control gorup at follow-up
## 
##  Two Sample t-test
## 
## data:  fup12_wmc_sum by tvc
## t = 7.813, df = 76, p-value = 2.508e-11
## alternative hypothesis: true difference in means between group 1 and group 2 is not equal to 0
## 95 percent confidence interval:
##  19.09369 32.15883
## sample estimates:
## mean in group 1 mean in group 2 
##       34.848485        9.222222
## 
##  Effect size (cohen's D)
## [1] 1.790619

HHS: Comparing mean HHS in both groups at baseline and follow-up

Baseline: The t-test was conducted to compare the baseline Harris hip scores of Group 1 and Group 2. The calculated t-value was 4.298 with 112 degrees of freedom. The associated p-value was <0.001, which is less than the chosen significance level of 0.05. Since the p-value is less than the chosen significance level (α = 0.05), we reject the null hypothesis. This indicates that there was significant difference in means of baseline Harris hip score between Group 1 and Group 2, with mean score being higher in group 1.The effect size for the difference was large (Cohen’s d > 0.8).

Follow-up: The t-test was conducted to compare the follow-up Harris hip scores of Group 1 and Group 2. The calculated t-value was -6.125 with 81 degrees of freedom. The associated p-value was <0.001, which is less than the chosen significance level of 0.05. Since the p-value is less than the chosen significance level (α = 0.05), we reject the null hypothesis. This indicates that there was significant difference in means of follow-up Harris hip score between Group 1 and Group 2, with mean score being higher in group 2. The effect size for the difference was large (Cohen’s d > 0.8).

## HHS in treatment and control gorup at baseline
## 
##  Two Sample t-test
## 
## data:  bl_hhs_total by tvc
## t = 4.2989, df = 112, p-value = 3.68e-05
## alternative hypothesis: true difference in means between group 1 and group 2 is not equal to 0
## 95 percent confidence interval:
##   6.934456 18.791898
## sample estimates:
## mean in group 1 mean in group 2 
##        68.20603        55.34286
## 
##  Effect size (cohen's D)
## [1] 0.8053741
## HHS in treatment and control group at follow-up
## 
##  Two Sample t-test
## 
## data:  fup_hhs_total by tvc
## t = -6.1251, df = 81, p-value = 3.097e-08
## alternative hypothesis: true difference in means between group 1 and group 2 is not equal to 0
## 95 percent confidence interval:
##  -19.64193 -10.00983
## sample estimates:
## mean in group 1 mean in group 2 
##        76.69079        91.51667
## 
##  Effect size (cohen's D)
## [1] 1.349445

ANCOVA

ANCOVA: Testing the difference in mean follow-up WOMAC scores

ANCOVA was carried out for testing the difference in mean follow-up WOMAC scores between treatment and control groups, while controlling/*adjusting for baseline differences in mean WOMAC scores

The independent variable tvc (Treatment) demonstrated a highly significant effect on the dependent variable (follow-up WOMAC score)(estimate = -24.7320, SE = 3.4066, t value = -7.260, p < 0.001***) i.e there was a difference of approximately -24.7 points in WOMAC scores of both groups after controlling for the baseline values.

The co-variate “baseline womac score” was marginally significant, suggesting a potential influence on the dependent variable.

# ANCOVA



model_ancova_wmc = lm (fup12_wmc_sum ~ bl_wmc_sum + tvc ,
              data = wp1_data) # Is the model selection okay?

summary(model_ancova_wmc)
## 
## Call:
## lm(formula = fup12_wmc_sum ~ bl_wmc_sum + tvc, data = wp1_data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -35.378  -7.369  -1.608   9.177  35.117 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  25.2901     5.2216   4.843 7.88e-06 ***
## bl_wmc_sum    0.2135     0.1233   1.731   0.0881 .  
## tvc2        -24.7320     3.4066  -7.260 5.18e-10 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 14.02 on 67 degrees of freedom
##   (50 observations deleted due to missingness)
## Multiple R-squared:  0.4445, Adjusted R-squared:  0.4279 
## F-statistic:  26.8 on 2 and 67 DF,  p-value: 2.803e-09
# Checking assumptions of the model
 
hist(residuals(model_ancova_wmc),
     col="darkgray")

## A histogram of residuals from a linear model.

plot(fitted(model_ancova_wmc),
     residuals(model_ancova_wmc))

ANCOVA:Testing the difference in mean follow-up Harris hip scores

ANCOVA was carried out for testing the difference in mean follow-up Harris hip scores between treatment and control groups, while controlling/*adjusting for baseline differences in mean Harris hip scores

The independent variable tvc (Treatment) demonstrated a highly significant effect on the dependent variable (follow-up Harris hip score)(estimate = 16.26183, SE = 2.72, t value = 5.96, p < 0.001***) i.e there was a difference of approximately 16 points in HHS scores of both groups after controlling for the baseline values.

The co-variate “baseline Harris hip sccore” was not significant, suggesting no influence on the dependent variable.

model_ancova_hhs = lm (fup_hhs_total ~ bl_hhs_total + tvc ,
              data = wp1_data)

summary(model_ancova_hhs)
## 
## Call:
## lm(formula = fup_hhs_total ~ bl_hhs_total + tvc, data = wp1_data)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -28.209  -7.508   2.061   5.047  21.244 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  72.17612    6.09499  11.842  < 2e-16 ***
## bl_hhs_total  0.05777    0.08192   0.705    0.483    
## tvc2         16.26183    2.72771   5.962  7.6e-08 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 10.99 on 75 degrees of freedom
##   (42 observations deleted due to missingness)
## Multiple R-squared:  0.3425, Adjusted R-squared:  0.3249 
## F-statistic: 19.53 on 2 and 75 DF,  p-value: 1.487e-07
# Checking assumptions of the model
 
hist(residuals(model_ancova_hhs),
     col="darkgray")

## A histogram of residuals from a linear model.

plot(fitted(model_ancova_hhs),
     residuals(model_ancova_hhs))

MIXED MODELS

Linear mixed model for testing differences in follow-up WOMAC scores

The mixed-effects model was fitted to investigate the impact of the treatment variable (tvc) on the outcome variable (WOMAC scores), accounting for potential interactions with time_point and individual-specific variations (id).

The results revealed that there was no difference in WOMAC scores between both groups at baseline (tvc2, Est. -1.30, p-value 0.64). The WOMAC score was significantly lower in group 2 at both 6 month (tvc2:time_pointfup6_wmc_sum, Est. -23.43, p-value < 0.001) and 12 month (tvc2:time_pointfup12_wmc_sum, Est. -26.27, p-value <0.001) follow up.

Pairwise comparisons revealed no significant differences in baseline WOMAC scores and scores at 6-month and 12-month follow-ups within Group 1. In Group 2, WOMAC scores were significantly lower at both the 6-month and 12-month follow-up assessments compared to the baseline scores within the same group. No significant differences were observed in baseline WOMAC scores between the two groups. However, at both the 6-month and 12-month follow-up time points, WOMAC scores were significantly lower in Group 2 compared to Group 1

#mixed model

l2 <- lmer(outcome ~ tvc * time_point + (1|id), data = long_data_wmc, REML= FALSE)
summ(l2)
Observations 262
Dependent variable outcome
Type Mixed effects linear regression
AIC 2093.14
BIC 2121.69
Pseudo-R² (fixed effects) 0.48
Pseudo-R² (total) 0.74
Fixed Effects
Est. S.E. t val. d.f. p
(Intercept) 41.11 1.99 20.69 194.85 0.00
tvc2 -1.30 2.76 -0.47 197.99 0.64
time_pointfup12_wmc_sum -4.46 2.37 -1.88 174.08 0.06
time_pointfup6_wmc_sum -4.41 2.20 -2.00 167.11 0.05
tvc2:time_pointfup12_wmc_sum -26.27 3.17 -8.28 171.36 0.00
tvc2:time_pointfup6_wmc_sum -23.43 3.07 -7.63 168.14 0.00
p values calculated using Satterthwaite d.f.
Random Effects
Group Parameter Std. Dev.
id (Intercept) 9.87
Residual 9.97
Grouping Variables
Group # groups ICC
id 110 0.49 
lsmeans(l2, pairwise~tvc * time_point)
## Registered S3 methods overwritten by 'broom':
##   method            from  
##   tidy.glht         jtools
##   tidy.summary.glht jtools
## $lsmeans
##  tvc time_point    lsmean   SE  df lower.CL upper.CL
##  1   bl_wmc_sum     41.11 2.01 201    37.15     45.1
##  2   bl_wmc_sum     39.81 1.94 208    35.99     43.6
##  1   fup12_wmc_sum  36.65 2.34 245    32.04     41.3
##  2   fup12_wmc_sum   9.07 2.05 224     5.03     13.1
##  1   fup6_wmc_sum   36.70 2.18 224    32.41     41.0
##  2   fup6_wmc_sum   11.97 2.09 229     7.86     16.1
## 
## Degrees-of-freedom method: kenward-roger 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                                estimate   SE  df t.ratio p.value
##  tvc1 bl_wmc_sum - tvc2 bl_wmc_sum         1.3012 2.79 204   0.466  0.9972
##  tvc1 bl_wmc_sum - tvc1 fup12_wmc_sum      4.4628 2.41 181   1.855  0.4332
##  tvc1 bl_wmc_sum - tvc2 fup12_wmc_sum     32.0369 2.87 213  11.157  <.0001
##  tvc1 bl_wmc_sum - tvc1 fup6_wmc_sum       4.4067 2.23 174   1.977  0.3598
##  tvc1 bl_wmc_sum - tvc2 fup6_wmc_sum      29.1400 2.90 216  10.062  <.0001
##  tvc2 bl_wmc_sum - tvc1 fup12_wmc_sum      3.1616 3.04 231   1.041  0.9037
##  tvc2 bl_wmc_sum - tvc2 fup12_wmc_sum     30.7356 2.14 175  14.376  <.0001
##  tvc2 bl_wmc_sum - tvc1 fup6_wmc_sum       3.1055 2.92 217   1.065  0.8947
##  tvc2 bl_wmc_sum - tvc2 fup6_wmc_sum      27.8388 2.17 176  12.808  <.0001
##  tvc1 fup12_wmc_sum - tvc2 fup12_wmc_sum  27.5740 3.11 236   8.858  <.0001
##  tvc1 fup12_wmc_sum - tvc1 fup6_wmc_sum   -0.0561 2.43 166  -0.023  1.0000
##  tvc1 fup12_wmc_sum - tvc2 fup6_wmc_sum   24.6772 3.14 238   7.870  <.0001
##  tvc2 fup12_wmc_sum - tvc1 fup6_wmc_sum  -27.6301 2.99 224  -9.230  <.0001
##  tvc2 fup12_wmc_sum - tvc2 fup6_wmc_sum   -2.8968 2.20 166  -1.317  0.7755
##  tvc1 fup6_wmc_sum - tvc2 fup6_wmc_sum    24.7333 3.02 227   8.197  <.0001
## 
## Degrees-of-freedom method: kenward-roger 
## P value adjustment: tukey method for comparing a family of 6 estimates

Linear mixed model for testing differences in follow-up Harris hip scores

The mixed-effects model was fitted to investigate the impact of the treatment variable (tvc) on the outcome variable (Harris hip scores), accounting for potential interactions with time_point and individual-specific variations (id). The results revealed that the baseline Harris hip scores were significantly lower in group2 when compared to group 1 (tvc2, Estimate = -12.90, SE = 2.61, t value = -4.93, p = <0.001). At follow-up, the Harris Hip score was significantly higher in group 2 when compared to group 1 (tvc2:time_pointfup_hhs_total , Estimate = 27.79, SE = 3.90, t value = 7.14, p = <0.001)

Pairwise comparisons revealed significant difference in baseline and follow-up HHS in group 1. The difference however was relatively smaller when compared to the difference in baseline and follow-up HHS in group 2.

l3 <- lmer(outcome ~ tvc * time_point + (1|id), data = long_data_hhs, REML= FALSE)
summ(l3)
Observations 197
Dependent variable outcome
Type Mixed effects linear regression
AIC 1609.05
BIC 1628.75
Pseudo-R² (fixed effects) 0.47
Pseudo-R² (total) 0.52
Fixed Effects
Est. S.E. t val. d.f. p
(Intercept) 68.23 1.83 37.23 195.70 0.00
tvc2 -12.90 2.61 -4.93 195.86 0.00
time_pointfup_hhs_total 8.23 2.81 2.92 98.37 0.00
tvc2:time_pointfup_hhs_total 27.79 3.90 7.14 96.20 0.00
p values calculated using Satterthwaite d.f.
Random Effects
Group Parameter Std. Dev.
id (Intercept) 4.07
Residual 13.35
Grouping Variables
Group # groups ICC
id 119 0.09 
lsmeans(l3, pairwise~tvc * time_point)
## $lsmeans
##  tvc time_point    lsmean   SE  df lower.CL upper.CL
##  1   bl_hhs_total    68.2 1.85 200     64.6     71.9
##  2   bl_hhs_total    55.3 1.89 200     51.6     59.0
##  1   fup_hhs_total   76.5 2.29 201     71.9     81.0
##  2   fup_hhs_total   91.3 2.11 201     87.2     95.5
## 
## Degrees-of-freedom method: kenward-roger 
## Confidence level used: 0.95 
## 
## $contrasts
##  contrast                                estimate   SE  df t.ratio p.value
##  tvc1 bl_hhs_total - tvc2 bl_hhs_total      12.90 2.64 200   4.882  <.0001
##  tvc1 bl_hhs_total - tvc1 fup_hhs_total     -8.23 2.85 112  -2.884  0.0240
##  tvc1 bl_hhs_total - tvc2 fup_hhs_total    -23.12 2.80 201  -8.242  <.0001
##  tvc2 bl_hhs_total - tvc1 fup_hhs_total    -21.13 2.97 201  -7.113  <.0001
##  tvc2 bl_hhs_total - tvc2 fup_hhs_total    -36.02 2.73 107 -13.199  <.0001
##  tvc1 fup_hhs_total - tvc2 fup_hhs_total   -14.89 3.12 201  -4.781  <.0001
## 
## Degrees-of-freedom method: kenward-roger 
## P value adjustment: tukey method for comparing a family of 4 estimates