Нахождение пикового значения в сигнале кривой колоколообразной формы с использованием R - PullRequest
1 голос
/ 18 марта 2020

У меня есть такие данные. Это для одного участника рулевого управления в течение 2,5 с. Меня интересует сигнал скорости рыскания (diff_YR)

df <- structure(list(ppid_trialn = c("1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39", 
"1_5_39", "1_5_39", "1_5_39", "1_5_39", "1_5_39"), diff_YR = c(0, 
0, 0, 0, 0, 0, 0, 0, 0, -0.0029370662937063, -0.00220279636363637, 
-0.0015501090909091, 0.00775058834498834, 0.00864802237762236, 
0.0260256428904429, 0.0375291375291375, 0.0518065249417249, 0.0751398564102564, 
0.0939044233100233, 0.108508151048951, 0.134778553939394, 0.149545457902098, 
0.164393939393939, 0.169452214452214, 0.180058275058275, 0.166678321678322, 
0.167494172494173, 0.169207459207459, 0.162272727272727, 0.157051282051282, 
0.158846153846154, 0.157377622377622, 0.14979020979021, 0.14513986013986, 
0.131025641025641, 0.112016317016317, 0.1104662004662, 0.0927622377622378, 
0.0775874125874126, 0.0607808857808858, 0.0456060606060607, 0.0353263403263403, 
0.0183566433566434, 0.00962703962703965, 0.00122377622377622, 
-0.00228438228438226, 0.00367132867132869, -0.00367132867132869, 
-0.0035897435897436, 0.000244755244755215, -0.00832167832167834, 
-0.0115034965034965, -0.012972027972028, -0.0214568764568764, 
-0.0257808857808857, -0.0317365967365967, -0.0406293706293706, 
-0.0457692307692308, -0.0470745920745922, -0.0545804195804196, 
-0.0703263403263403, -0.0793822843822844, -0.0917832167832168, 
-0.109568764568765, -0.117808857808858, -0.131025641025641, -0.143834498834499, 
-0.146363636363636, -0.150198135198135, -0.164801873193473, -0.171573416083916, 
-0.179242416083916, -0.186993001165501, -0.187400928904429, -0.177529136363636, 
-0.171410259254079, -0.171981353379953, -0.160151519347319, -0.153706300699301, 
-0.147587422377622, -0.152156168531469, -0.153379953613054, -0.158682984382284, 
-0.16357808974359, -0.157703964335664, -0.161701633799534, -0.174755243076923, 
-0.171655011655012, -0.178100233100233, -0.17997668997669, -0.190501165501166, 
-0.195559440559441, -0.206002331002331, -0.223135198135198, -0.215874125874126, 
-0.229825174825175, -0.23472027972028, -0.230559440559441, -0.232762237762238, 
-0.225174825174825, -0.212365967365967, -0.196783216783217, -0.18462703962704, 
-0.168881118881119, -0.132983682983683, -0.109160839160839, -0.0757925407925408, 
-0.0495221445221446, -0.0265151515151514, 0.00155011655011648, 
0.0209673659673659, 0.0438111888111889, 0.0527039627039627, 0.0650233100233101, 
0.079055944055944, 0.0835431235431236, 0.0921911421911422, 0.10467365967366, 
0.110547785547786, 0.119277389277389, 0.128414918414918, 0.14, 
0.148648018648019, 0.161456876456876, 0.167167832167832, 0.182424242424242, 
0.201841491841492), timestamp_zero = c(0.300099999999986, 0.316699999999997, 
0.33359999999999, 0.350200000000001, 0.36669999999998, 0.383299999999991, 
0.400000000000006, 0.416599999999988, 0.433300000000003, 0.449999999999989, 
0.466800000000006, 0.483299999999986, 0.5, 0.516699999999986, 
0.5334, 0.550099999999986, 0.566699999999997, 0.583399999999983, 
0.599999999999994, 0.616600000000005, 0.633299999999991, 0.650099999999981, 
0.666599999999988, 0.683300000000003, 0.699999999999989, 0.7166, 
0.733299999999986, 0.75, 0.766699999999986, 0.7834, 0.800099999999986, 
0.816699999999997, 0.833499999999987, 0.849999999999994, 0.86669999999998, 
0.883399999999995, 0.900000000000006, 0.916699999999992, 0.933300000000003, 
0.949999999999989, 0.966700000000003, 0.983299999999986, 1, 1.01669999999999, 
1.0334, 1.04999999999998, 1.0667, 1.08329999999998, 1.1001, 1.11660000000001, 
1.13329999999999, 1.15000000000001, 1.16659999999999, 1.1833, 
1.19999999999999, 1.2166, 1.23319999999998, 1.25, 1.26669999999999, 
1.2833, 1.29999999999998, 1.31659999999999, 1.33329999999998, 
1.34999999999999, 1.36669999999998, 1.38339999999999, 1.40000000000001, 
1.41669999999999, 1.4333, 1.45009999999999, 1.4667, 1.48329999999999, 
1.5, 1.51669999999999, 1.5334, 1.54999999999998, 1.56659999999999, 
1.58329999999998, 1.59999999999999, 1.61660000000001, 1.63339999999999, 
1.65000000000001, 1.66659999999999, 1.6833, 1.69999999999999, 
1.7166, 1.73329999999999, 1.7501, 1.76659999999998, 1.7834, 1.79999999999998, 
1.81659999999999, 1.83329999999998, 1.8501, 1.86660000000001, 
1.88319999999999, 1.8999, 1.91659999999999, 1.9333, 1.94989999999999, 
1.9667, 1.98339999999999, 1.9999, 2.01659999999998, 2.03319999999999, 
2.04990000000001, 2.06659999999999, 2.08320000000001, 2.09989999999999, 
2.11660000000001, 2.13319999999999, 2.15000000000001, 2.16649999999998, 
2.1832, 2.19979999999998, 2.2165, 2.23310000000001, 2.2499, 2.26659999999998, 
2.28319999999999, 2.29999999999998, 2.31649999999999, 2.33320000000001, 
2.34989999999999, 2.3665, 2.38319999999999, 2.3998), SWA = c(0, 
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0.0900001, 0.1800001, 
0.36, 0.45, 0.81, 0.9900001, 1.35, 1.8, 2.16, 2.61, 3.15, 3.51, 
3.96, 4.41, 4.77, 5.13, 5.67, 6.03, 6.39, 6.75, 7.11, 7.47, 7.65, 
7.92, 8.1, 8.28, 8.46, 8.46, 8.55, 8.55, 8.55, 8.55, 8.55, 8.55, 
8.55, 8.55, 8.55, 8.46, 8.46, 8.46, 8.28, 8.28, 8.19, 8.01, 7.92, 
7.74, 7.56, 7.47, 7.2, 6.84, 6.57, 6.21, 5.85, 5.49, 5.13, 4.68, 
4.23, 3.78, 3.33, 2.79, 2.34, 1.89, 1.44, 1.08, 0.72, 0.27, -0.1799999, 
-0.45, -0.9, -1.26, -1.8, -2.16, -2.61, -3.06, -3.42, -3.96, 
-4.5, -4.95, -5.49, -6.03, -6.57, -7.2, -7.83, -8.46, -8.91, 
-9.63, -10.08, -10.62, -11.16, -11.52, -11.88, -12.15, -12.42, 
-12.51, -12.51, -12.51, -12.51, -12.42, -12.24, -12.06, -11.88, 
-11.61, -11.43, -11.16, -10.89, -10.62, -10.17, -9.81, -9.54, 
-9.09, -8.64, -8.19, -7.65)), class = c("tbl_df", "tbl", "data.frame"
), row.names = c(NA, -127L))

При построении этого сигнала во времени получается следующий график:

ggplot(df, aes(x = timestamp_zero, y = diff_YR)) +
  geom_line()

plot

Я хочу сделать, это выбрать значения (то есть diff_YR, время, угол поворота рулевого колеса (SWA) и т. Д. c) на пике сигнала diff_YR. Сначала я использовал следующий код:

upperthreshold = 0.13 # upper threshold for consistent steering response
lowerthreshold = 0.05 # lower threshold for when response initiated

steering_peak <- df %>%
  dplyr::select(ppid_trialn, diff_YR, timestamp_zero, SWA) %>%
  group_by(ppid_trialn) %>%
  dplyr::mutate(peakYaw = max(diff_YR)) %>%
  dplyr::filter(max(diff_YR) > upperthreshold, min(diff_YR) < lowerthreshold) %>%
  slice(1:max(which(diff_YR > upperthreshold, 1))) %>%
  slice(max(which(diff_YR == peakYaw))) %>%
  ungroup() %>%
  transmute(ppid_trialn, PeakYawRateChange = diff_YR, PeakSteeringTime = timestamp_zero, peakSWA = SWA)

Вышеупомянутый код группируется пробным способом. Затем он определяет максимальный сигнал diff_YR. Затем я фильтрую сигнал, чтобы включить наибольшее значение diff_YR через верхний порог и самое маленькое значение diff_YR ниже нижнего порога. Затем я разрезаю вектор diff_YR от первого значения до максимального значения diff_YR, которое выше верхнего порога. Затем я разрезаю строку данных, когда diff_YR равен максимальному рысканию.

Это прекрасно работает, если мое максимальное значение diff_YR, однако, как вы можете видеть из графика, это не дело для этого испытания. Поэтому я в конечном итоге выбираю «пик», который не является истинным пиком сигнала. Одно из предложений, которое у меня было, было отфильтровать временной интервал, где происходит пик. Однако для некоторых испытаний пик наступает позднее во времени, поэтому мне нужно все это при рассмотрении нескольких испытаний.

Следовательно, мой вопрос: есть ли у кого-нибудь идея, как я могу определить и отфильтровать только колокольчик? отклик кривой для каждого из моих испытаний, а затем выберите строку данных в моем кадре данных, которая соответствует моему истинному пику сигнала diff_YR.

Большое спасибо заранее!

1 Ответ

0 голосов
/ 19 марта 2020

Вы можете попробовать использовать пакеты, которые находят пики и позволяют определить порог и т. Д. c, например, ниже, я использую findpeaks из pracma, где вы можете указать несколько параметров, таких как минимальная высота пика и минимальная пиковое расстояние. Я перебираю несколько настроек для minpeakpeakdistance, потому что проще установить что-то для минимальной высоты пика:

library(pracma)
plotPeaks = function(df,i){
p = findpeaks(df$diff_YR,minpeakdistance=i,minpeakheight=0)
plot(df$timestamp_zero,df$diff_YR,main=paste("minpeakdistance=",i))
abline(v=df$timestamp_zero[p[,2]])
}

par(mfrow=c(2,2))
for(I in c(3,5,7,9)){
plotPeaks(df,I)
}

enter image description here

Теперь, если вы объедините минимум На расстоянии пика 7 с минимальной высотой пика, скажем, 0,1, вы можете получить что-то разумное.

Приведенный выше пик поиска не учитывает время, и, как вы можете видеть, он также подвержен ударам. Вы также можете сначала сгладить данные или взять скользящее среднее, чтобы выровнять их:

df$mean_diff_YR = rollmean(df$diff_YR,5,fill=NA)
p = findpeaks(df$mean_diff_YR,minpeakdistance=7,minpeakheight=0)
plot(df$timestamp_zero,df$diff_YR,main=paste("minpeakdistance=",7))
abline(v=df$timestamp_zero[p[,2]])

enter image description here

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...