The packages that we are going to use
library(ggplot2)
library(dplyr)
df<-read.table("rawdata_27Oct.csv", header=TRUE, sep='\t')
See the structure of the data
str(df)
## 'data.frame': 564 obs. of 7 variables:
## $ ID : int 2 3 6 9 10 17 21 24 25 26 ...
## $ Group : Factor w/ 2 levels "A","B": 1 1 1 1 1 1 1 1 1 1 ...
## $ Score : int 10 10 6 8 8 8 7 9 9 7 ...
## $ Type : Factor w/ 2 levels "NRS","ODI 50": 1 1 1 1 1 1 1 1 1 1 ...
## $ Phase : Factor w/ 5 levels "12-Month","3-Month",..: 4 4 4 4 4 4 4 4 4 4 ...
## $ Crossover: Factor w/ 3 levels "","NO","YES": 2 2 2 2 2 2 2 2 2 2 ...
## $ Dropout : Factor w/ 3 levels "","NO","YES": 2 2 2 3 3 3 2 2 2 2 ...
Change the Order of the Phase Levels
df$Phase<-factor(df$Phase, levels=c( "Baseline" , "Month", "3-Month", "6-Month", "12-Month" ))
Analysis for the ODI-50 Type
model<-lm(Score~Group+Phase, data=subset(df, Type=='ODI 50'))
summary(model)
##
## Call:
## lm(formula = Score ~ Group + Phase, data = subset(df, Type ==
## "ODI 50"))
##
## Residuals:
## Min 1Q Median 3Q Max
## -20.4899 -6.0208 0.7479 5.9792 18.9468
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 22.783 1.284 17.738 < 2e-16 ***
## GroupB 4.238 1.114 3.803 0.000183 ***
## PhaseMonth -6.531 1.616 -4.042 7.23e-05 ***
## Phase3-Month -5.795 1.631 -3.552 0.000463 ***
## Phase6-Month -5.590 1.684 -3.319 0.001051 **
## Phase12-Month -7.968 2.054 -3.879 0.000137 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 8.549 on 231 degrees of freedom
## (45 observations deleted due to missingness)
## Multiple R-squared: 0.1488, Adjusted R-squared: 0.1304
## F-statistic: 8.077 on 5 and 231 DF, p-value: 4.85e-07
ANOVA for the ODI-50 Type
anova(model)
## Analysis of Variance Table
##
## Response: Score
## Df Sum Sq Mean Sq F value Pr(>F)
## Group 1 1142.8 1142.85 15.6364 0.0001021 ***
## Phase 4 1808.9 452.22 6.1872 9.556e-05 ***
## Residuals 231 16883.6 73.09
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Testing also the interaction between Groups and Phase
model<-lm(Score~Group+Phase+Group*Phase, data=subset(df, Type=='ODI 50'))
anova(model)
## Analysis of Variance Table
##
## Response: Score
## Df Sum Sq Mean Sq F value Pr(>F)
## Group 1 1142.8 1142.85 15.7193 9.848e-05 ***
## Phase 4 1808.9 452.22 6.2200 9.125e-05 ***
## Group:Phase 4 379.9 94.97 1.3063 0.2685
## Residuals 227 16503.7 72.70
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Analysis for the NRS Type
model<-lm(Score~Group+Phase, data=subset(df, Type=='NRS'))
summary(model)
##
## Call:
## lm(formula = Score ~ Group + Phase, data = subset(df, Type ==
## "NRS"))
##
## Residuals:
## Min 1Q Median 3Q Max
## -6.5349 -1.4830 0.1766 1.7291 5.7291
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 7.4830 0.3766 19.872 < 2e-16 ***
## GroupB 1.3700 0.3261 4.201 3.79e-05 ***
## PhaseMonth -3.2121 0.4740 -6.776 1.01e-10 ***
## Phase3-Month -2.6596 0.4786 -5.557 7.50e-08 ***
## Phase6-Month -2.3181 0.4941 -4.691 4.64e-06 ***
## Phase12-Month -2.8076 0.5941 -4.726 3.98e-06 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 2.508 on 232 degrees of freedom
## (44 observations deleted due to missingness)
## Multiple R-squared: 0.2429, Adjusted R-squared: 0.2265
## F-statistic: 14.88 on 5 and 232 DF, p-value: 1.149e-12
ANOVA for the NRS Type
anova(model)
## Analysis of Variance Table
##
## Response: Score
## Df Sum Sq Mean Sq F value Pr(>F)
## Group 1 116.84 116.843 18.578 2.411e-05 ***
## Phase 4 351.17 87.793 13.959 3.214e-10 ***
## Residuals 232 1459.15 6.289
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Testing also the interaction between Groups and Phase
model<-lm(Score~Group+Phase+Group*Phase, data=subset(df, Type=='NRS'))
anova(model)
## Analysis of Variance Table
##
## Response: Score
## Df Sum Sq Mean Sq F value Pr(>F)
## Group 1 116.84 116.843 19.0520 1.929e-05 ***
## Phase 4 351.17 87.793 14.3152 1.928e-10 ***
## Group:Phase 4 60.86 15.214 2.4808 0.04475 *
## Residuals 228 1398.29 6.133
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Test if at the beggining the two samples were equivalent
A_NRS<-subset(df, Type=="NRS"&Phase=="Baseline"&Group=="A", select=Score )
B_NRS<-subset(df, Type=="NRS"&Phase=="Baseline"&Group=="B", select=Score )
A_ODI<-subset(df, Type=="ODI 50"&Phase=="Baseline"&Group=="A", select=Score )
B_ODI<-subset(df, Type=="ODI 50"&Phase=="Baseline"&Group=="B", select=Score )
t.test(A_NRS, B_NRS)
##
## Welch Two Sample t-test
##
## data: A_NRS and B_NRS
## t = 0.015626, df = 52.704, p-value = 0.9876
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -0.6321410 0.6420665
## sample estimates:
## mean of x mean of y
## 8.230769 8.225806
t.test(A_ODI, B_ODI)
##
## Welch Two Sample t-test
##
## data: A_ODI and B_ODI
## t = -0.55542, df = 48.41, p-value = 0.5812
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -4.991784 2.830494
## sample estimates:
## mean of x mean of y
## 24.50000 25.58065
Test if in a Month the two samples were equivalent
A_NRS<-subset(df, Type=="NRS"&Phase=="Month"&Group=="A", select=Score )
B_NRS<-subset(df, Type=="NRS"&Phase=="Month"&Group=="B", select=Score )
A_ODI<-subset(df, Type=="ODI 50"&Phase=="Month"&Group=="A", select=Score )
B_ODI<-subset(df, Type=="ODI 50"&Phase=="Month"&Group=="B", select=Score )
t.test(A_NRS, B_NRS)
##
## Welch Two Sample t-test
##
## data: A_NRS and B_NRS
## t = -2.3939, df = 52.029, p-value = 0.02031
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.4313256 -0.3020078
## sample estimates:
## mean of x mean of y
## 4.000000 5.866667
t.test(A_ODI, B_ODI)
##
## Welch Two Sample t-test
##
## data: A_ODI and B_ODI
## t = -2.8226, df = 52.995, p-value = 0.006694
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -11.426797 -1.933203
## sample estimates:
## mean of x mean of y
## 14.92 21.60
Test if in 3 Months the two samples were equivalent
A_NRS<-subset(df, Type=="NRS"&Phase=="3-Month"&Group=="A", select=Score )
B_NRS<-subset(df, Type=="NRS"&Phase=="3-Month"&Group=="B", select=Score )
A_ODI<-subset(df, Type=="ODI 50"&Phase=="3-Month"&Group=="A", select=Score )
B_ODI<-subset(df, Type=="ODI 50"&Phase=="3-Month"&Group=="B", select=Score )
t.test(A_NRS, B_NRS)
##
## Welch Two Sample t-test
##
## data: A_NRS and B_NRS
## t = -2.1508, df = 50.738, p-value = 0.03628
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -3.0273400 -0.1040886
## sample estimates:
## mean of x mean of y
## 4.720000 6.285714
t.test(A_ODI, B_ODI)
##
## Welch Two Sample t-test
##
## data: A_ODI and B_ODI
## t = -1.8308, df = 51, p-value = 0.07297
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -8.8354413 0.4068699
## sample estimates:
## mean of x mean of y
## 17.00000 21.21429
Test if in 6 Months the two samples were equivalent
A_NRS<-subset(df, Type=="NRS"&Phase=="6-Month"&Group=="A", select=Score )
B_NRS<-subset(df, Type=="NRS"&Phase=="6-Month"&Group=="B", select=Score )
A_ODI<-subset(df, Type=="ODI 50"&Phase=="6-Month"&Group=="A", select=Score )
B_ODI<-subset(df, Type=="ODI 50"&Phase=="6-Month"&Group=="B", select=Score )
t.test(A_NRS, B_NRS)
##
## Welch Two Sample t-test
##
## data: A_NRS and B_NRS
## t = -3.6015, df = 43.249, p-value = 0.0008107
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -4.353443 -1.228376
## sample estimates:
## mean of x mean of y
## 4.409091 7.200000
t.test(A_ODI, B_ODI)
##
## Welch Two Sample t-test
##
## data: A_ODI and B_ODI
## t = -2.7643, df = 44.866, p-value = 0.008248
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -12.088217 -1.897238
## sample estimates:
## mean of x mean of y
## 15.72727 22.72000
Test if in 12 Months the two samples were equivalent
A_NRS<-subset(df, Type=="NRS"&Phase=="12-Month"&Group=="A", select=Score )
B_NRS<-subset(df, Type=="NRS"&Phase=="12-Month"&Group=="B", select=Score )
A_ODI<-subset(df, Type=="ODI 50"&Phase=="12-Month"&Group=="A", select=Score )
B_ODI<-subset(df, Type=="ODI 50"&Phase=="12-Month"&Group=="B", select=Score )
t.test(A_NRS, B_NRS)
##
## Welch Two Sample t-test
##
## data: A_NRS and B_NRS
## t = -0.31203, df = 23.697, p-value = 0.7577
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -2.721014 2.006728
## sample estimates:
## mean of x mean of y
## 5.142857 5.500000
t.test(A_ODI, B_ODI)
##
## Welch Two Sample t-test
##
## data: A_ODI and B_ODI
## t = -0.20964, df = 18.787, p-value = 0.8362
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## -9.849463 8.057255
## sample estimates:
## mean of x mean of y
## 16.28571 17.18182
A table of the Average and Standard Deviation By Group, Period and Type
df%>%group_by(Group, Phase, Type)%>%summarise(Average=mean(Score, na.rm=TRUE), Standard_Dev=sd(Score, na.rm=TRUE))
## Source: local data frame [20 x 5]
## Groups: Group, Phase [?]
##
## Group Phase Type Average Standard_Dev
## (fctr) (fctr) (fctr) (dbl) (dbl)
## 1 A Baseline NRS 8.230769 1.210213
## 2 A Baseline ODI 50 24.500000 7.925907
## 3 A Month NRS 4.000000 2.813657
## 4 A Month ODI 50 14.920000 7.921069
## 5 A 3-Month NRS 4.720000 2.590367
## 6 A 3-Month ODI 50 17.000000 7.895146
## 7 A 6-Month NRS 4.409091 2.737032
## 8 A 6-Month ODI 50 15.727273 8.339132
## 9 A 12-Month NRS 5.142857 2.983471
## 10 A 12-Month ODI 50 16.285714 9.193882
## 11 B Baseline NRS 8.225806 1.175009
## 12 B Baseline ODI 50 25.580645 6.515490
## 13 B Month NRS 5.866667 2.956388
## 14 B Month ODI 50 21.600000 9.629695
## 15 B 3-Month NRS 6.285714 2.706058
## 16 B 3-Month ODI 50 21.214286 8.862560
## 17 B 6-Month NRS 7.200000 2.549510
## 18 B 6-Month ODI 50 22.720000 8.997778
## 19 B 12-Month NRS 5.500000 2.844452
## 20 B 12-Month ODI 50 17.181818 11.600157
Conclusions
- At the beginning of the experiment both groups were equivalent which means that the samples are comparable. We tested both pain scale tpyes
- There was statistically significant difference between Group A and B when the treatment started
- In all phases the score was statistical significant different than from the start of the experiment (i.e. Baseline)
- In a Month the group A was significantly better than group B in both NRS and ODI-50 at 5% level of significance (P-Values are 0.02031 and 0.006694 respectively)
- In 3 Months the group A was significantly better than group B in both NRS and ODI-50 at 5% level of significance (P-Values are 0.03628 and 0.07297 respectively)
- In 6 Months the group A was significantly better than group B in both NRS and ODI-50 at 1% level of significance (P-Values are 0.001 and 0.008 respectively)
- In 12 Months both groups are equivalent (P-Values are 0.7577 and 0.8362 respectively)
Some relevant Boxplot charts
qplot(Group, Score, data=subset(df, Type=='NRS'), geom = "boxplot", color=Phase, main="NRS by Group and Phase")
qplot(Group, Score, data=subset(df, Type=='ODI 50'), geom = "boxplot", color=Phase, main="ODI 50 by Group and Phase")
qplot(Group, Score, data=df, geom = "boxplot", color=Group, facets=Type~Phase, main="Group and Phase per Type")
qplot(Group, Score, data=subset(df, Type=='NRS'), geom = "boxplot", color=Group, facets=.~Phase, main="NRS by Group and Phase")
qplot(Group, Score, data=subset(df, Type=='ODI 50'), geom = "boxplot", color=Group, facets=.~Phase, main="ODI 50 by Group and Phase")
Plot of the Average Pain Score of NRS and ODI-50
ggplot(subset(df, Type=='NRS'), aes(x=Group, y=Score)) + stat_summary(fun.y="mean", geom="bar")+facet_grid(.~Phase)+aes(color=Group)+ggtitle("NRS Average Score")
ggplot(subset(df, Type=='ODI 50'), aes(x=Group, y=Score)) + stat_summary(fun.y="mean", geom="bar")+facet_grid(.~Phase)+aes(color=Group)+ggtitle("ODI 50 Average Score")
ggplot(df, aes(x=Group, y=Score)) + stat_summary(fun.y="mean", geom="bar")+facet_grid(Type~Phase)+aes(color=Group)+ggtitle("Average Score By Period and Type")
