БПФ со сглаживанием окна Хэмминга и извлечением данных из автоматически созданного графика - PullRequest
Новая
       123

БПФ со сглаживанием окна Хэмминга и извлечением данных из автоматически созданного графика

0 голосов
/ 14 июля 2020

Мне нужна помощь в анализе временных рядов, в частности, быстрые преобразования Фурье со сглаживанием окна Хэмминга.

TL; DR

  1. is fftper () подходящую функцию для БПФ со сглаживанием окна Хэмминга в R?

  2. Как я могу извлечь или сгенерировать значения частоты выходного графика fftper?

  3. Я ищу суточные циклы, поэтому, чтобы преобразовать данные частоты обратно в переменную «время», нужно ли мне разделить значения частоты на 1/24? (например, здесь пример кори: http://web.stanford.edu/class/earthsys214/notes/series.html)

Длинная версия

У меня есть временная метка acousti c обнаружение данные для группы отдельных животных. Для каждого человека я отсортировал количество обнаружений по часам и преобразовал их во временной ряд с помощью ts () с частотой 24 (просматривая дневные шаблоны в обнаружениях). С этими данными я хочу применить сглаживание окна Хэмминга, а затем БПФ и создать периодограмму этих данных. Я также хочу извлечь значения частоты (ось x) и преобразовать их из частоты в период времени (часы).

Мне удалось сгенерировать периодограмму на БПФ и обработать оконные данные с помощью функции fftper() в формате TSSS package. Автоматически созданные графики выглядят правильно. Теперь я хочу извлечь значения частоты (значения оси x) и значения мощности (значения оси y), используемые в графике, чтобы я мог преобразовать данные частоты (значения оси x) обратно во временную переменную (т.е. используя, как мне кажется, частота / (1/24)?), А затем красиво построить его с помощью ggplot. fftper генерирует объект spg, который имеет такую ​​структуру: 

List of 4
 $ period         : num [1:65] 10.43 1.95 2.36 2.57 1.9 ...
 $ smoothed.period: num [1:65] 0.815 0.601 0.364 0.374 0.348 ...
 $ log.scale      : chr "TRUE"
 $ tsname         : chr "hourly_ts"
 - attr(*, "class")= chr "spg"

Я могу извлечь значения оси Y (сглаженные значения периода или мощности) с помощью FFTpower <- FFT[["smoothed.period"]], но я не вижу, где значения оси x сохранены или выясните, как их сгенерировать.

Есть идеи? Заранее спасибо!

Фиктивные данные: 

#Data
df <- read.table(text =
               "timestampUTC    ID
'2017-10-02 19:23:27'    47280
'2017-10-02 19:26:48'    47280
'2017-10-02 19:27:23'    47280
'2017-10-02 19:31:46'    47280
'2017-10-02 23:52:15'    47280
'2017-10-02 23:53:26'    47280
'2017-10-02 23:55:13'    47280
'2017-10-03 19:53:50'    47280
'2017-10-03 19:55:23'    47280
'2017-10-03 19:58:26'    47280
'2017-10-04 13:15:13'    47280
'2017-10-04 13:16:42'    47280
'2017-10-04 13:21:39'    47280
'2017-10-04 19:34:54'    47280
'2017-10-04 19:55:28'    47280
'2017-10-04 20:08:23'    47280
'2017-10-04 20:21:43'    47280
'2017-10-05 04:55:48'    47280
'2017-10-05 04:57:04'    47280
'2017-10-05 05:18:40'    47280
'2017-10-07 21:24:19'    47280
'2017-10-07 21:25:36'    47280
'2017-10-07 21:29:25'    47280", header = T)

Код: 

#convert datetime
df$timestampUTC<-as.POSIXct(df$timestampUTC, format = "%Y-%m-%d %H:%M:%S", tz="UTC")

#keep only datetime column and add second column with frequency of 1
df<-df %>%
  select(timestampUTC)
df<-data.frame(df,Frequency=1)

#bin into hours
hourly_detections <- df %>%
  mutate(processed_hour = floor_date(timestampUTC, "hour")) %>%
  group_by(processed_hour)%>%
  summarise(count = sum(Frequency))

#set time frame using max and min hours
time_frame <- as_datetime(c(min(floor_date(df$timestampUTC,"hour")),(max(ceiling_date(df$timestampUTC,"hour"))-1)),tz="Australia/Sydney")

#combine detection hour and non-detections hour dfs
all_hours <- data.frame(processed_hour = seq(time_frame[1], time_frame[2], by = "hour"))

#build df with every hour and set count to 0 for 'new' hours
hourly_detections <- hourly_detections %>%
  right_join(all_hours, by = "processed_hour") %>%
  mutate(count = ifelse(test = is.na(count),yes  = 0,no   = count))
hourly_detections<-hourly_detections[order(hourly_detections$processed_hour),]

#set up time series
hourly_ts <- ts(hourly_detections$count, start= min(hourly_detections$processed_hour), frequency=24) 

#FFT with hamming widow smoothing
FFT<-fftper(hourly_ts, window = 2, plot = TRUE)

#extract y (power) values
FFTPower<-FFT[["smoothed.period"]]

#extract x values?
...