The Systolic Blood Pressure Dataset

data(sbp)
head(sbp)

##   meth item repl   y
## 1    J    1    1 100
## 2    J    2    1 108
## 3    J    3    1  76
## 4    J    4    1 108
## 5    J    5    1 124
## 6    J    6    1 122

#sbp <- sbp %>% rename("item" = Id, "repl" = Rep, "meth" = Obs, "y" = Vic) %>% dplyr::select(-Sub)
sbp <- sbp %>% mutate(item = factor(item)) %>% filter(meth != "R")

sbp.wide <-(sbp %>% spread(meth,y)) %>% group_by(item) %>% summarize(J = mean(J,na.rm=TRUE),S=mean(S,na.rm=TRUE))
sbp.mean <-(sbp %>% spread(meth,y)) %>% group_by(item) %>% summarize(J = mean(J,na.rm=TRUE),S=mean(S,na.rm=TRUE))

mod.1 <- lm( y ~ meth , data = sbp)
summary(mod.1)

confint(mod.1)

mod.2 <- lm( S ~ J , data = (sbp %>% spread(meth,y)))
summary(mod.2)

confint(mod.2)

mod.3 <- lm( S ~ J , data = (sbp.mean))
summary(mod.3)

confint(mod.3)

ttest(mod.3, 2,1)

mod.4 <- lme( y ~ meth,random = ~ 1|item, data = sbp)

summary(mod.4)


# tidy(mod.4)

mod.4 <- lmer( y ~ meth + (1|item), data = sbp)

summary(mod.4)


tidy(mod.4)

mod.5 <- lme( S ~ J , random = ~ 1|item, data = sbp.wide)

summary(mod.5)


#intervals(mod.5)

mod.5 <- lme( S ~ J , random = ~ 1|item, data = (sbp %>% spread(meth,y)))

summary(mod.5)

intervals(mod.5)

mod.5 <- lmer( S ~ J + (1|item), data = (sbp %>% spread(meth,y)))

summary(mod.5)

sbp.2 <- sbp %>% spread(meth,y)
mod.5 <- lmer( S ~ J + (1|item), data = sbp.2)

#summary(mod.5)

tidy(mod.5)

#mod.6 <-  lmer( y ~ meth-1 + (1|item) + (1|meth) , data = sbp,
                control = lmerControl(
                           optimizer ='optimx', optCtrl=list(method='nlminb')))
#summary(mod.6)

#tidy(mod.6)

sbp.aug <- augment(mod.6,sbp)
head(sbp.aug)

var.test(.wtres~ meth, data = sbp.aug)

We conclude that there is no significant difference in within-item error terms.

grouping.levels(mod.6, "item")
mod.6.inf  <- influence(mod.6,group="item")

mod.6.inf

CDs<- cooks.distance(mod.6.inf)
subjectnames<- rownames(cooks.distance(mod.6.inf))
CDdf <- data.frame(subjectnames,CDs)

dfbeta <- dfbetas(mod.6.inf)

subjectnames <- rownames(dfbetas(mod.6.inf))

DFBetadf <- data.frame(subjectnames,dfbeta)

head(DFBetadf)

library(ggplot2)
ggplot(data=DFBetadf,aes(x=methJ,y=methS)) + geom_point() + geom_abline() + 
  stat_smooth(method="lm",colour="red",se=FALSE)