In this example we will show how you can use the autohrf package to automatically generate data-informed event models for GLM modeling of task-based fMRI data. Let us start this example by loading required libraries and the data from the spatial working memory experiment (SWM).
# libs
library(autohrf)
# load the data
df <- swm
head(df)## roi t y
## 1 L_2_ROI 0 2.445844
## 2 L_46_ROI 0 2.155406
## 3 L_6a_ROI 0 4.436729
## 4 L_7Am_ROI 0 1.169667
## 5 L_AIP_ROI 0 4.132354
## 6 L_FEF_ROI 0 7.725104The loaded data frame has 504 observations, each with 3 variables (roi, t, and y) roi denotes the region of interest, t the time stamp and y the value of the BOLD signal. Note that input data for the autohrf package should be always organized in this manner.
Next, we will define two different model constraints which the autohrf function will use to find the best fitting model given these constraints. Note that this is only an illustrative example and the set amount of the population size and iterations of the genetic algorithm is way too small for any kind of proper analysis.
# 2 events: delay, response
model2 <- data.frame(event = c("delay", "response"),
start_time = c(0, 10),
end_time = c(10, 15))
# 3 events: encoding, delay, response
model3 <- data.frame(event = c("encoding", "delay", "response"),
start_time = c(0, 0, 10),
end_time = c(0.5, 10, 15))Once we define the model constraints we can use the our models we can use the autohrf function to automatically find model parameters that fit our data best.
model_constraints <- list(model2, model3)
# to speed vignette building we here load results from a previous autohrf run
autofit <- swm_autofit
# in practice you should run
#autofit <- autohrf(df, model_constraints, tr = 2.5)When the automated fitting process is completed, we can use the get_best_models function to extract the best model for each of the provided constraints.
best <- get_best_models(autofit)##
## ----------------------------------------
##
## Model 1
##
## Fitness: 0.8893632
##
## event start_time duration
## 1 delay 0.30000 3.844288
## 2 response 10.93357 4.066429
##
## ----------------------------------------
##
## ----------------------------------------
##
## Model 2
##
## Fitness: 0.9347359
##
## event start_time duration
## 1 encoding 0.39 0.110000
## 2 delay 0.46 9.540000
## 3 response 11.88 3.085662
##
## ----------------------------------------Based on calculated fitness scores we can see that the second model fits our data better. Furthermore, we can use the plot_best_models function to visually analyze the best models.
plot_best_models(autofit)