计算时间序列的速度(时间的一阶导数)

Calculate speed (first derivative in time) of time series

我想根据 xts 时间序列中的数据计算速度。我的数据看起来像这样(下面的命令生成样本,实际数据要大得多)。

measurement <- xts(c(7.9, 8.6, 12.7, 13.3), 
                   as.POSIXct(c("2016-02-24 06:00:00",
                                "2016-02-24 07:30:00",
                                "2016-02-24 09:15:00",
                                "2016-02-24 11:15:00")))
names(measurement) <- "pos"

head(measurement)
                     pos
2016-02-24 06:00:00  7.9
2016-02-24 07:30:00  8.6
2016-02-24 09:15:00 12.7
2016-02-24 11:15:00 13.3

如果我使用 diff,我会得到连续值之间的变化。

diff(measurement)
                    pos
2016-02-24 06:00:00  NA
2016-02-24 07:30:00 0.7
2016-02-24 09:15:00 4.1
2016-02-24 11:15:00 0.6

但是,我想考虑时间和计算速度(例如时间段的变化)。我尝试过使用以下肮脏的道路。

measurement$time  <- as.numeric(index(measurement))
measurement$speed  <- diff(measurement$pos) / diff(measurement$time) * 3600

head(measurement)
                    pos  time       speed
2016-02-24 06:00:00  7.9 1456290000        NA
2016-02-24 07:30:00  8.6 1456295400 0.4666667
2016-02-24 09:15:00 12.7 1456301700 2.3428571
2016-02-24 11:15:00 13.3 1456308900 0.3000000

一定有更简单、更优雅的方法来做到这一点。请记住,我对 R 还是个新手,可能会遗漏一些东西。

我不知道 base R 或 xts 包提供了一种直接计算物理速率的方法。在数小时内进行 speed 计算的更直接的方法可能是

# alternative time conversion 
  del_t <- diff(index(measurement))
  units(del_t) <- "hours"
  measurement$speed_hr <- diff(measurement$pos, na.pad=FALSE)/as.numeric(del_t)

可以指定 return 由 diff 编辑的时间单位,因此在这种情况下可以设置为 hours

更一般地说,您的速度计算在数据点处存在不连续性。在某些情况下,还需要 speed 跨数据点的连续性。 speed 然后可以使用 R 的 splinefun 例程计算,该例程不仅可以 return 位置的样条插值,还可以给出一阶导数,给出速度的近似值,它既连续又取决于超过两个相邻的数据点。代码可能看起来像

# provides continuity of speed at data points
# perform cubic spline fit
  spline_fit <- splinefun(x=index(measurement), y=measurement$pos, method="natural")

# add spline speeds to measurements
  measurement$spline_spd <- spline_fit(index(measurement), deriv=1)*3600

样条速度与原始计算的速度不一致,但这似乎是数据点连续性受限的结果。一个情节可能有助于澄清这一点。

# make a sequence of plot times for spline fits
  spline_times <- seq(start(measurement), end(measurement), length.out=40)
# plot positions
  par(mfcol=c(2,1))
  plot(spline_times, spline_fit(spline_times, deriv=0), col="red", type="b", 
       main="Positions", xlab = "Time", ylab="Pos")
  points(index(measurement), measurement$pos, type="p", pch=19, cex=1.1)
  legend("topleft", legend = c("-- pos data","-- spline pos interpolations"), 
         text.col = c("black","red"), y.intersp=.2, yjust=3., bty="n", cex=1.3)
# plot speeds
  plot(spline_times, spline_fit(spline_times, deriv=1)*3600, col="red", type="b",
       main="Speeds", xlab="Time", ylab="Speed")
  lines(index(measurement), measurement$speed, type="b", pch=19, cex=1.1)
  legend("topleft", legend = c("-- speed calculation","-- spline speed interpolations"), 
       text.col = c("black","red"), y.intersp=.2, yjust=3., bty="n", cex=1.3)

四个数据点对于良好的样条拟合来说确实太少了,但也许这传达了一般的想法。