Analysis of “Fyfe NH dogs_OU_whole day avgs_STATS_2020.xlsx” pre-post trend.
First read the data and make transformations
setwd("C:/Users/marti/OneDrive/Documents/andras/Kristin")
library(tidyverse)
library(readxl)
library(lme4)
library(lsmeans)
library(ggplot2)
dts<-read_xlsx("Fyfe NH dogs_OU_whole day avgs_STATS_2020.xlsx",sheet = "All_NH_dogs")
dts<-mutate(dts,aft= (1+(`days pre and post injection`<=0)),pre=c("post","pre")[aft])
Let’s look at trends pre-post injection
ggplot(dts,aes(x=`days pre and post injection`,y=`Whole day averages`))+
geom_point()+geom_smooth(method = "lm")+
facet_grid(rows = vars(Name),cols = vars(pre))
## `geom_smooth()` using formula 'y ~ x'
It does seem suspicious, let’s break data into two parts, keep zero in both subsets and re-do plots separately
g1<-filter(dts,`days pre and post injection`<=0)
g2<-filter(dts,`days pre and post injection`>=0)
ggplot(g1,aes(x=`days pre and post injection`,y=`Whole day averages`))+
geom_point()+geom_smooth(method = "lm")+
facet_wrap(~Name)+ggtitle("Pre")
## `geom_smooth()` using formula 'y ~ x'
ggplot(g2,aes(x=`days pre and post injection`,y=`Whole day averages`))+
geom_point()+geom_smooth(method = "lm")+
facet_wrap(~Name)+ggtitle("Post")
## `geom_smooth()` using formula 'y ~ x'
Some dogs have large variation between eyes. let’s confirm sample size and do averaging.
table(dts$Name,dts$`days pre and post injection`)
##
## -9 -6 0 5 33 34 37
## Claremont 2 2 2 2 2 2 2
## Merrimack 2 2 2 2 2 2 2
## Nashua 2 2 2 2 2 2 2
## Sunny 2 2 2 2 2 0 0
## Winnie 0 0 2 2 2 0 0
confirmed: two obs per dog, but different timepoings. Look at averaged data now.
avg<-group_by(dts,Name,`days pre and post injection`)
dtm<-summarize(avg,y=mean(`Whole day averages`),pre=c("post","pre")[mean(aft)])
## `summarise()` has grouped output by 'Name'. You can override using the `.groups` argument.
ggplot(dtm,aes(x=`days pre and post injection`,y=y))+
geom_point()+geom_smooth(method = "lm")+
facet_grid(rows = vars(Name),cols = vars(pre))
## `geom_smooth()` using formula 'y ~ x'
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in qt((1 - level)/2, df): NaNs produced
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
## Warning in max(ids, na.rm = TRUE): no non-missing arguments to max; returning -
## Inf
g1<-filter(dtm,`days pre and post injection`<=0)
g2<-filter(dtm,`days pre and post injection`>=0)
m1<-lmer(y~`days pre and post injection`+(1|Name),data=g1)
summary(m1)
## Linear mixed model fit by REML ['lmerMod']
## Formula: y ~ `days pre and post injection` + (1 | Name)
## Data: g1
##
## REML criterion at convergence: 56.8
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -1.74953 -0.38665 -0.01138 0.38196 1.57579
##
## Random effects:
## Groups Name Variance Std.Dev.
## Name (Intercept) 0.6666 0.8165
## Residual 4.4792 2.1164
## Number of obs: 13, groups: Name, 5
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 16.1968 0.9886 16.383
## `days pre and post injection` -0.0326 0.1545 -0.211
##
## Correlation of Fixed Effects:
## (Intr)
## `papinjctn` 0.707
m2<-lmer(y~`days pre and post injection`+(1|Name)+(`days pre and post injection`-1|Name),data=g2)
summary(m2)
## Linear mixed model fit by REML ['lmerMod']
## Formula:
## y ~ `days pre and post injection` + (1 | Name) + (`days pre and post injection` -
## 1 | Name)
## Data: g2
##
## REML criterion at convergence: 117.6
##
## Scaled residuals:
## Min 1Q Median 3Q Max
## -2.2371 -0.5039 0.2396 0.5666 1.1683
##
## Random effects:
## Groups Name Variance Std.Dev.
## Name (Intercept) 0.6209 0.7880
## Name.1 `days pre and post injection` 0.3006 0.5482
## Residual 5.4253 2.3292
## Number of obs: 21, groups: Name, 5
##
## Fixed effects:
## Estimate Std. Error t value
## (Intercept) 15.1937 0.8684 17.497
## `days pre and post injection` 0.1308 0.2476 0.528
##
## Correlation of Fixed Effects:
## (Intr)
## `papinjctn` -0.093