Last active
July 21, 2022 20:10
-
-
Save WillForan/4469423e225c75e0ed1c8a3c5813f5ac to your computer and use it in GitHub Desktop.
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
library(ggplot) | |
library(gamm4) | |
library(dplyr) | |
library(checkmate) # for expect_subset dataframe names test | |
met_res_for_roi <- function(MRSI_input, this_roi, met_name, minage=18, maxage=25) { | |
# have columns we need | |
expect_subset(c("roi", "age", "dateNumeric", met_name), names(MRSI_input)) | |
roi_met <- MRSI_input %>% | |
filter(roi %in% this_roi, age >=minage, age <=maxage) %>% | |
mutate(met = .data[[met_name]], # map whatever metabolite we're using (met_name) onto the 'met' column | |
gamresids = gam(met ~ s(age, k=3), data=.) %>% residuals, | |
dateAdjust = gam(gamresids ~ s(dateNumeric, k=3), data=.) %>% predict, | |
met_adj = met - dateAdjust) | |
return(roi_met) | |
} | |
plot_met_adjusted <- function(roi_met) { | |
# have columns we need | |
expect_subset(c("age", "gamresids", "dateAdjust", "met", "met_adj", "visitnum", "subjID"), | |
names(roi_met)) | |
qassert(roi_met$met_adj, "n+") # nonzero numeric vector in met_adj | |
cowplot::plot_grid( | |
ggplot(roi_met) + aes(x=age, y=gamresids) + geom_smooth(method="loess") + geom_point(), | |
ggplot(roi_met) + aes(x=vdate, y=dateAdjust) + geom_point() + geom_smooth(method="loess"), | |
ggplot(roi_met) + aes(x=vdate, y=met_adj, color=visitnum, group=visitnum) + geom_point() + geom_smooth(method="loess"), | |
ggplot(roi_met) + aes(x=age, y=met_adj) + geom_point() + stat_smooth(method='loess', se=F) + geom_path(aes(group = subjID)), | |
ggplot(roi_met) + aes(x=age, y=met) + geom_point() + stat_smooth(method='loess', se=F) + geom_path(aes(group = subjID)) | |
) | |
} | |
test_data <- data.frame(subjID=1:32, roi=1, age=rep(10:25,2), gaba=randu[1:32,1], dateNumeric=sample(rep(1:8,4)),visitnum=1) %>% mutate(vdate=dateNumeric) | |
adj_df <- met_res_for_roi(test_data, 1, "gaba") | |
plot_met_adjusted(adj_df) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment