#Q1)What assumption must we test to include a variable as a blocking factor? 
#A1)To include a variable as a blocking factor, we must test the assumption of homogeneity of variance. This means that we need to check if the variance of the dependent variable is equal across all levels of the blocking factor. If this assumption is not met, then it may indicate that the blocking factor is affecting the variance of the dependent variable and therefore needs to be included in the analysis to control for its effect.


#Q2)Recognize the IV, DV, block and create a table for the following research statement - “A company is planning to investigate the motor skills of elderly population. The company separates the target population into three age categories: 60 – 69, 70 – 79, and above 80 then randomly assign the participants in the study to one of the three task conditions. After individuals have completed the task, their performance will be compared.”

#A2)Independent Variable (IV): Age category (60-69, 70-79, and above 80)
#Dependent Variable (DV): Performance on the motor skill task
#Block: Conditions: 1, 2 ,3 


#Q3) Use the data “Lab 3” with the research question to perform a fine report. 
# H0: Motor skills (performance score) of elderly population are equal 
#H1: At least one pair of scores differ from one another. 
library("readxl")
library("dplyr")
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
library("moments")
library("pastecs")
## 
## Attaching package: 'pastecs'
## The following objects are masked from 'package:dplyr':
## 
##     first, last
library("compute.es")

motor = read_excel("C:/Users/apoor/Downloads/Lab3.xlsx")
motor = motor %>% mutate(Age = Age %>% as.factor(), Condition = Condition %>% as.factor())

#Checking overall distribution 
plot(density(motor$Performance_score))

qqnorm(motor$Performance_score)

#D'Agostino skeness test

agostino.test(motor$Performance_score)
## 
##  D'Agostino skewness test
## 
## data:  motor$Performance_score
## skew = -0.11171, z = -0.45976, p-value = 0.6457
## alternative hypothesis: data have a skewness
plot(density(motor$Performance_score))

lm = lm(Performance_score ~ Age, data=motor) 
res = resid(lm)
plot(motor$Performance_score, res, ylab="Residuals", xlab="Performance")
abline(0, 0)

#Based on the residual plot, the residuals do not look normally distributed, but are independent.

bartlett.test(motor$Performance_score, motor$Age)
## 
##  Bartlett test of homogeneity of variances
## 
## data:  motor$Performance_score and motor$Age
## Bartlett's K-squared = 1.0587, df = 2, p-value = 0.589
tapply(motor$Performance_score, motor$Age, var)
##        1        2        3 
## 12.89901 18.99570 15.83744

#We can see here that the variances are closeby accross the age groups.

m <- aov(Performance_score ~ Condition*Age  , data =motor)
summary(m)
##               Df Sum Sq Mean Sq F value Pr(>F)    
## Condition      2 1199.0   599.5 313.667 <2e-16 ***
## Age            2 1549.6   774.8 405.389 <2e-16 ***
## Condition:Age  4   22.6     5.7   2.961 0.0246 *  
## Residuals     80  152.9     1.9                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
m2 <- aov(Performance_score ~ Age+Condition  , data = motor)
summary(m2)
##             Df Sum Sq Mean Sq F value Pr(>F)    
## Age          2 1549.7   774.9   370.8 <2e-16 ***
## Condition    2 1198.9   599.5   286.9 <2e-16 ***
## Residuals   84  175.5     2.1                   
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

#We can see that there is no significant interaction between the treatment and the block.

pairwise.t.test(motor$Performance_score, motor$Age, paired = FALSE)
## 
##  Pairwise comparisons using t tests with pooled SD 
## 
## data:  motor$Performance_score and motor$Age 
## 
##   1       2      
## 2 3.5e-05 -      
## 3 3.2e-15 4.8e-07
## 
## P value adjustment method: holm
TukeyHSD(m2)
##   Tukey multiple comparisons of means
##     95% family-wise confidence level
## 
## Fit: aov(formula = Performance_score ~ Age + Condition, data = motor)
## 
## $Age
##           diff        lwr       upr p adj
## 2-1  -4.512792  -5.403864 -3.621720     0
## 3-1 -10.310345 -11.216147 -9.404543     0
## 3-2  -5.797553  -6.688625 -4.906480     0
## 
## $Condition
##          diff       lwr       upr p adj
## 2-1 -4.189467 -5.095808 -3.283126     0
## 3-1 -9.005432 -9.904688 -8.106176     0
## 3-2 -4.815965 -5.699331 -3.932599     0
t = by(motor$Performance_score, motor$Age, stat.desc)
mes(t$`1`['mean'], t$`2`['mean'], t$`1`['std.dev'], t$`2`['std.dev'], t$`1`['nbr.val'], t$`2`['nbr.val'], verbose = FALSE)$d
## [1] 1.13
mes(t$`1`['mean'], t$`3`['mean'], t$`1`['std.dev'], t$`3`['std.dev'], t$`1`['nbr.val'], t$`3`['nbr.val'], verbose = FALSE)$d
## [1] 2.72

#Summary: We used a one way ANOVA. Since all assumptions were met, no adjustments were made. The results suggested that the age has a significant effect on the performance {F (2, 89) = 370.8, p < .001}. A Bonferroni hoc test shows that there exists significant difference between the three groups. We then did a Tukey’s hoc test and the result suggested that there is a significant difference between the groups. The effect was large, difference of score between 60 – 69 and 70 – 79 is Cohen’ D = 1.13. and difference of score between 60 – 69 and above 80 is Cohen’ D = 2.72. So we reject the H0.