-
Notifications
You must be signed in to change notification settings - Fork 1
/
Copy pathhrv_stats.Rmd
125 lines (98 loc) · 4.03 KB
/
hrv_stats.Rmd
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
---
title: "HRV Statistics"
author: "Daniel N. Albohn"
date: "11/08/2017"
output:
html_document:
theme: cosmo
css: style.css
editor_options:
chunk_output_type: console
---
```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE, warning = FALSE, message = FALSE, error = FALSE, eval = FALSE)
```
We finally have our cleaned IBI time series (derived from the cleaned and annotated ECG time series)!
From here, we can use the `RHRV` package to analyze our data in a number of ways. Below, I build
the time series and perform a frequency analysis on it. Finally, we will extract the high frequency
signal for two of the marked points of interest.
First, we load the data and create a new, blank HRV class using `CreateHRVData()`.
This is where we will store everything as we proceed.
```{r}
library(RHRV)
# Set up some global variables
path <- 'data/hrv_tutorial/'
file <- 'sub1101_ecg_clean.txt'
name <- sub("*_ecg_clean.txt", "", file)
hrv.data = CreateHRVData()
hrv.data = SetVerbose(hrv.data, FALSE)
hrv.data = LoadBeatRR(hrv.data, RecordName=file.path(path,file), RecordPath=".", scale = .001)
```
Below, we load the trigger file we created in the previous step
and extract the info we need to overlay on the time series so the
program knows where we want it to derive the statistics from. We add
the triggers using `AddEpisodes()`.
```{r}
# We add the info about the episodes
file_ev <- sub("*_ecg_clean.txt", "", file)
load(file.path(path,paste0(name,"_trigger.RData")))
InitTime <- episodes$InitTime
Type <- episodes$trigger
Duration <- episodes$Duration
Value <- episodes$Value
hrv.data = AddEpisodes(hrv.data, InitTimes = episodes$InitTime,
Tags = episodes$Type,
Durations = episodes$Duration,
Values = episodes$Value)
```
Next, we can derive the instantaneous heart rate time series by using the cleaned IBI
time series. This is possible because "instantaneous heart rate can be defined as the
inverse of the time separation between two consecutive heart beats."
```{r}
hrv.data = BuildNIHR(hrv.data)
hrv.data = FilterNIHR(hrv.data)
# plot all tags
# png(filename = paste("data/plots/",name,"_tagged_plot.png",sep=""), width=1000, height=669,
# units="px")
PlotNIHR(hrv.data, Tag=episodes$Type)
# dev.off()
hrv.data = InterpolateNIHR(hrv.data, freqhr = 4)
```

Once everything is loaded into our `hrv.data` object, we can perform
the frequency analysis with two commands: one to build a
frequency analysis sub-object in our `hrv.data` object, and another
to actual calculate the powerband. We can plot the powerband to check for
any anomalies.
```{r}
#Perform frequency analysis
##Calculating spectrogram and power per band using wavelet analysis:
hrv.data = CreateFreqAnalysis(hrv.data)
hrv.data = CalculatePowerBand(hrv.data, indexFreqAnalysis = 1, type="wavelet",
wavelet="d4", bandtolerance=0.1)
# plot powerband for all files
# png(filename = paste("data/plots/",name,"_powerband.png",sep=""), width=1000, height=669,
# units="px")
PlotPowerBand(hrv.data, normalized = TRUE, hr = TRUE, Tag = "all")
# dev.off()
```

Once we are satisfied with the frequency analysis, we can begin to chop up
the data into parts that we will use to compare.
```{r}
# Save the data by stimulus type:
splitting.data1 = SplitPowerBandByEpisodes(hrv.data, indexFreqAnalysis = 1, Tag = c("event_type_1.29"))
Baseline <- log(mean(splitting.data1$OutEpisodes$HF))
splitting.data2 = SplitPowerBandByEpisodes(hrv.data, indexFreqAnalysis = 1, Tag = c("event_type_1.92"))
Task <- log(mean(splitting.data2$OutEpisodes$HF))
subject_nr <- readr::parse_number(file)
sub <- cbind.data.frame(subject_nr,Baseline,Task)
write.table(sub, file = "data/hrv_tutorial/hrv_extract_data.csv", sep = ",", append = FALSE,
col.names = TRUE, row.names = FALSE)
# Clean up
rm(list = ls())
```
Let's make sure that our data were saved properly.
```{r eval=TRUE}
head(read.csv('data/hrv_tutorial/hrv_extract_data.csv'))
```