library(PupillometryR)
#> Loading required package: 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
#> Loading required package: ggplot2
data("pupil_data")
#Check that IDs are not numeric
pupil_data$ID <- as.character(pupil_data$ID)
#remove participant number 8, who had problematic data
pupil_data <- subset(pupil_data, ID != 8)
#blinks were registered as -1, so replace with NAs
pupil_data$LPupil[pupil_data$LPupil == -1] <- NA
pupil_data$RPupil[pupil_data$RPupil == -1] <- NA
Sdata <- make_pupillometryr_data(data = pupil_data,
subject = ID,
trial = Trial,
time = Time,
condition = Type)
new_data <- replace_missing_data(data = Sdata,
rpupil = RPupil,
lpupil = LPupil)
#> Warning in replace_missing_data(data = Sdata, rpupil = RPupil, lpupil
#> = LPupil): replace_missing_data will only help if you have missing
#> timepoints, and a reliable time column.
plot(new_data, pupil = LPupil, group = 'condition')
#> Warning: Removed 3639 rows containing non-finite values (stat_summary).
regressed_data <- regress_data(data = new_data,
pupil1 = RPupil,
pupil2 = LPupil)
mean_data <- calculate_mean_pupil_size(data = regressed_data, rpupil = RPupil, lpupil = LPupil)
filtered_data <- filter_data(data = mean_data,
pupil = mpupil,
filter = 'hanning',
degree = 11)
int_data <- interpolate_data(data = filtered_data,
pupil = mpupil,
type = 'linear')
#plot
plot(int_data, pupil = mpupil, group = 'subject')
base_data <- baseline_data(data = int_data,
pupil = mpupil,
start = 0,
stop = 100)
plot(base_data, pupil = mpupil, group = 'subject')
library(mgcv)
#> Loading required package: nlme
#>
#> Attaching package: 'nlme'
#> The following object is masked from 'package:dplyr':
#>
#> collapse
#> This is mgcv 1.8-26. For overview type 'help("mgcv-package")'.
base_data$IDn <- as.numeric(base_data$ID)
base_data$Typen <- ifelse(base_data$Type == 'Easy', .5, -.5)
base_data$Trialn <- as.numeric(substr(base_data$Trial, 5, 5))
base_data$Trialn <- ifelse(base_data$Typen == .5, base_data$Trialn, base_data$Trialn + 3)
base_data$ID <- as.factor(base_data$ID)
model_data <- base_data
m1 <- bam(mpupil ~ s(Time) + s(Time, by = Typen),
data = base_data,
family = gaussian)
library(itsadug)
#> Loading required package: plotfunctions
#>
#> Attaching package: 'plotfunctions'
#> The following object is masked from 'package:ggplot2':
#>
#> alpha
#> Loaded package itsadug 2.3 (see 'help("itsadug")' ).
diffs <- get_difference(
m1,
comp = list(Typen = c(-0.5, 0.5)),
cond = list(Time = unique(base_data$Time)),
print.summary = TRUE)
#> Summary:
#> * Time : numeric predictor; with 600 values ranging from 16.666667 to 10000.000000.
ggplot(diffs) +
aes(x = Time, y = difference) +
geom_hline(yintercept = 0, color = "white", size = 2) +
geom_ribbon(
aes(ymin = difference - CI, ymax = difference + CI),
alpha = .3) +
geom_line() +
ggtitle("This confidence interval is probably too narrow.") +
labs(subtitle = "The model should probably have by-subject random smooths.")
Created on 2019-02-26 by the reprex package (v0.2.1)