Skip to content

Commit

Permalink
u
Browse files Browse the repository at this point in the history
  • Loading branch information
@DominiqueMakowski
DominiqueMakowski committed Aug 20, 2024
1 parent 8a2443e commit c7ab9b5
Show file tree
Hide file tree
Showing 3 changed files with 150 additions and 33 deletions.
16 changes: 8 additions & 8 deletions analysis/1_physio.py
View file
Original file line number Diff line number Diff line change
@@ -1,8 +1,5 @@
import io
import os
import time

import autoreject
import matplotlib.pyplot as plt
import mne
import neurokit2 as nk
Expand All @@ -11,15 +8,17 @@
import PIL
import pyllusion as ill
import scipy.stats
import urllib.request
import requests

mne.set_log_level(verbose="WARNING")

# Convenience functions ======================================================================
r = urllib.request.urlopen(
"https://raw.githubusercontent.com/RealityBending/PrimalsInteroception/main/analysis/0_clean_physio.py"
).read()
eval(compile(r.decode(), "<string>", "single"))
# Download the load_physio() function
exec(
requests.get(
"https://raw.githubusercontent.com/RealityBending/PrimalsInteroception/main/analysis/func_load_physio.py"
).text
)


def qc_physio(df, info, sub, plot_ecg=[], plot_ppg=[]):
Expand Down Expand Up @@ -164,6 +163,7 @@ def qc_physio(df, info, sub, plot_ecg=[], plot_ppg=[]):
).statistic

# Hear Rate Variability (HRV) -------------------------------------------------------------
print(" - HCT - HRV")

hrv = pd.concat(
[
Expand Down
161 changes: 139 additions & 22 deletions analysis/2_cleaning.qmd
View file
Original file line number Diff line number Diff line change
Expand Up @@ -37,56 +37,173 @@ The initial sample consisted of `r report::report_participants(df, age="Age", ge
### MAIA

```{r}
#| code-fold: false
make_hist <- function(data, scales="fixed") {
data |>
pivot_longer(cols = everything()) |>
ggplot(aes(x=value, fill=name)) +
geom_histogram(bins=30) +
theme_minimal() +
scale_fill_viridis_d(guide="none") +
facet_wrap(~name, scales=scales)
}
select(df, starts_with("MAIA"), -matches("\\d$"), -contains("_R")) |>
make_hist()
```

select(df, starts_with("MAIA")) |>
select(-contains("_R"), -matches("\\d$")) |>
estimate_density() |>
ggplot(aes(x=x, y=y, color=Parameter)) +
geom_line()
```{r}
make_correlation <- function(data) {
cor <- data |>
correlation(method="pearson", redundant = TRUE, p_adjust="none") |>
correlation::cor_sort() |>
correlation::cor_lower()
cor |>
mutate(val = paste0(insight::format_value(r), format_p(p, stars_only=TRUE))) |>
mutate(Parameter2 = fct_rev(Parameter2)) |>
# mutate(Parameter1 = fct_relabel(Parameter1, \(x) str_remove_all(x, "Feedback_")),
# Parameter2 = fct_relabel(Parameter2, \(x) str_remove_all(x, "Feedback_"))) |>
ggplot(aes(x=Parameter1, y=Parameter2)) +
geom_tile(aes(fill = r), color = "white") +
geom_text(aes(label = val), size = 3) +
labs(title = "Correlation Matrix") +
scale_fill_gradient2(
low = "#2196F3",
mid = "white",
high = "#F44336",
breaks = c(-1, 0, 1),
guide = guide_colourbar(ticks=FALSE),
midpoint = 0,
na.value = "grey85",
limit = c(-1, 1)) +
theme_minimal() +
theme(legend.title = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.text.x = element_text(angle = 45, hjust = 1))
}
select(df, starts_with("MAIA"), -matches("\\d$"), -contains("_R")) |>
make_correlation()
```


### IAS

```{r}
#| code-fold: false
df$IAS_Total <- rowSums(select(df, starts_with("IAS"))) / 21
select(df, starts_with("IAS")) |>
select(-contains("_R"), -matches("\\d$")) |>
estimate_density() |>
ggplot(aes(x=x, y=y, color=Parameter)) +
geom_line()
select(df, starts_with("IAS"), -matches("\\d$"), -contains("_R")) |>
make_hist()
```

### HCT

```{r}
#| code-fold: false
select(df, starts_with("HCT")) |>
estimate_density() |>
ggplot(aes(x=x, y=y, color=Parameter)) +
geom_line()
make_hist()
```

```{r}
select(df, starts_with("HCT")) |>
make_correlation()
```

### PI-99

```{r}
#| code-fold: false
select(df, starts_with("PI")) |>
select(-matches("\\d")) |>
# Remove cols containing digit
estimate_density() |>
ggplot(aes(x=x, y=y, color=Parameter)) +
geom_line()
geom_line(linewidth=1) +
theme_minimal()
```


### HRV

```{r}
select(df, starts_with("HRV")) |>
make_hist(scales="free")
```


:::


## Exclusions


```{r}
#| code-fold: false
outliers <- list()
df$HRV_AI[df$HRV_AI > 52] <- NA
df$HRV_MeanNN[df$HRV_MeanNN > 1100] <- NA
df$HRV_MeanNN_HCT[df$HRV_MeanNN_HCT > 1100] <- NA
df$HRV_RMSSD[df$HRV_RMSSD > 100] <- NA
df$HRV_RMSSD_HCT[df$HRV_RMSSD_HCT > 100] <- NA
df$HRV_SDNN[df$HRV_SDNN > 150] <- NA
df$HRV_SDNN_HCT[df$HRV_SDNN_HCT > 150] <- NA
```


## Quick Analysis


```{r}
#| code-fold: false
make_cor2 <- function(data1, data2) {
cor <- correlation(data1, data2, method="pearson", p_adjust="none")
cor |>
mutate(val = paste0(insight::format_value(r), format_p(p, stars_only=TRUE))) |>
mutate(Parameter2 = fct_rev(Parameter2)) |>
# mutate(Parameter1 = fct_relabel(Parameter1, \(x) str_remove_all(x, "Feedback_")),
# Parameter2 = fct_relabel(Parameter2, \(x) str_remove_all(x, "Feedback_"))) |>
ggplot(aes(x=Parameter1, y=Parameter2)) +
geom_tile(aes(fill = r), color = "white") +
geom_text(aes(label = val), size = 3) +
labs(title = "Correlation Matrix") +
scale_fill_gradient2(
low = "#2196F3",
mid = "white",
high = "#F44336",
breaks = c(-1, 0, 1),
guide = guide_colourbar(ticks=FALSE),
midpoint = 0,
na.value = "grey85",
limit = c(-1, 1)) +
theme_minimal() +
theme(legend.title = element_blank(),
axis.title.x = element_blank(),
axis.title.y = element_blank(),
axis.text.x = element_text(angle = 45, hjust = 1))
}
make_cor2(
select(df, starts_with("HCT")),
select(df, starts_with("MAIA"), starts_with("IAS"), -matches("\\d$"), -contains("_R"))
)
make_cor2(
select(df, starts_with("HRV")),
select(df, starts_with("MAIA"), starts_with("IAS"), -matches("\\d$"), -contains("_R"))
)
make_cor2(
select(df, starts_with("HCT")),
select(df, starts_with("HRV"))
)
```

## Save

```{r}
Expand Down
6 changes: 3 additions & 3 deletions analysis/func_load_physio.py
View file
Original file line number Diff line number Diff line change
Expand Up @@ -29,7 +29,6 @@ def load_physio(path, sub):
# [10000, 10000 + rs.info["sfreq"] * 60 * 8], ymin=0, ymax=10, color="red"
# )
events = {"onset": [10000], "duration": [rs.info["sfreq"] * 60 * 8]}
rs.to_data_frame().plot(subplots=True)

assert len(events["onset"]) == 1 # Check that there is only one event

Expand Down Expand Up @@ -87,9 +86,10 @@ def consecutive_nans(ch="AF7"):

assert len(events["onset"]) == 6 # Check that there are 6 epochs (the 6 intervals)


# Interpolate signal interruptions
if sub in ["sub-13"]:
hct = hct.apply_function(nk.signal_fillmissing, picks="PPG_Muse", method="backward")
hct = hct.apply_function(
nk.signal_fillmissing, picks="PPG_Muse", method="backward"
)

return rs, hct

0 comments on commit c7ab9b5

Please sign in to comment.